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V 

This  report  describes  19  investigations  in  the  field  of  seismic 
discrimination.  The  contributions  are  grouped  as  follows:  6  de¬ 
scribe  progress  in  the  development  of  a  Seismic  Data  Center; 
3  describe  research  into  seismic  algorithms  that  will  be  imple¬ 
mented  at  the  Data  Center;  7  are  concerned  with  the  problem  of 
source  characterization  using  the  seismic  moment-tensor  for¬ 
mulation;  and  3  describe  recent  improvements  in  our  in-house 
computer  systems. 
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This  is  the  thirty-second  Semiannual  Technical  Summary  (SATS)  report  describing  the 
activities  of  Lincoln  Laboratory  funded  under  Project  Vela  Uniform.  This  report  covers  the 
period  1  October  1979  to  il  March  1980.  Project  Vela  Uniform  is  a  program  of  research  into 
the  discrimination  between  earthquakes  and  nuclear  explosions  by  seismic  means.  An  impor¬ 
tant  recent  emphasis  of  the  project  is  in  the  development  of  the  data-handling  and  analysis  tech¬ 
niques  that  may  be  appropriate  for  the  monitoring  of  a  potential  Comprehensive  Test  Ban  Treaty, 
presently  under  negotiation.  During  FY  1979,  Lincoln  Laboratory  completed  a  design  for  a 
Seismic  Data  Center  (SDC)  that  would  satisfy  this  requirement,  and  the  FY  1980  program  con¬ 
sists  of  beginning  the  actual  development  of  such  a  data  center,  and  carrying  out  seismic  re¬ 
search  with  special  emphasis  on  those  areas  directly  related  to  the  operations  of  the  SDC. 

The  Lincoln  Laboratory  program  is  focusing  on  the  development  of  three  subsystems  of  the 
SDC  -  namely  the  Seismic  Analysis  Station  (SAS),  the  Local  Computer  Network,  and  the  Data 
Base  Management  System.  Initial  work  has  focused  heavily  on  the  SAS  and,  in  particular,  on  a 
pre-prototype  software  system  referred  to  as  Mod  I.  This  system,  running  on  our  in-house 
PDP-ll/45  computer,  will  allow  experimentation  and  development  of  new-generation  techniques 
for  seismic -data  analysis  and  event  bulletin  preparation,  that  will  be  implemented  on  the  SAS 
hardware  later  in  the  current  fiscal  year. 

The  Data  Base  Management  System  and  the  Local  Computer  Network  are  in  the  design  and 
procurement  phases.  The  Data  Base  Management  System  design  has  been  refined  substantially 
from  that  described  in  the  SDC  Design  report,  and  now  includes  a  more  comprehensive  discussion 
of  the  indexing  of  waveform  data  for  rapid  retrieval.  The  Local  Computer  Network  requirements 
have  been  formulated  in  detail,  and  are  presently  being  procured. 

We  have  begun  to  examine  the  seismic  algorithms  that  will  be  employed  in  the  SDC.  Auto¬ 
matic  association  and  location  algorithms  have  been  selected  for  implementation  in  the  Mod  I 
system.  As  experience  in  the  use  of  this  system  grows,  we  expect  to  improve  or  replace  these 
algorithms  with  better  ones.  We  have  continued  our  interest  in  detection  algorithms,  and  de¬ 
scribe  an  algorithm  that  is  more  sensitive  to  emergent  arrivals.  A  synthetic  program  to  gen¬ 
erate  both  artificial  seismicity  and  arrivals  at  a  prescribed  set  of  stations  is  described.  This 
will  be  used  to  test  and  evaluate  SDC  algorithms. 

A  major  emphasis  of  the  seismic  research  program  has  been  in  the  area  of  source  character¬ 
ization  using  the  moment-tensor  formulation.  Attempts  have  been  made  to  understand  the  gener¬ 
ation  of  anomalous  Rayleigh  waves  by  certain  explosions  at  the  Eastern  Kazakh  test  site.  It  is 
shown  that  certain  components  of  the  moment  tensor  are  very  poorly  resolved  for  events  with 
zero  depth.  Analysis  of  the  anomalous  events  strongly  suggests  that  they  were  accompanied  by 
faulting.  A  new  method  for  the  simultaneous  inversion  of  surface -wave  data  for  both  moment - 
tensor  and  source  location  is  described.  This  method  shows  considerable  promise  for  the  anal¬ 
ysis  of  special  events. 

We  continue  to  develop  our  in-house  computer  systems.  New  enhancements  to  the  UNIX 
operating  system,  a  new  graphics  package,  and  the  installation  of  UNIX-compatible  Fortran  77 
are  described. 


M.A.  Chinnery 


SEISMIC  DISCRIMINATION 


I.  SEISMIC  DATA  CENTER  DEVELOPMENT 

\.  I>  U'hC.  ROUND  AND  CURRENT  STATUS 

On  10  September  1979,  Lincoln  Laboratory  completed  a  Special  Internal  Report  to 
DARPA  NMRO  entitled  ’'Design  of  a  Seismic  Data  Center."  This  design  study  has  been  reviewed 
and  i  opted  by  DA R PA/ NMRO  as  a  working  basis  for  the  development  of  a  Seismic  Data  Cen¬ 
ter  (S DC),  which  will  fulfill  the  following  objectives: 

(1)  Support  research  and  development  directed  at  developing  a  new  genera¬ 
tion  of  seismic -data  management  and  analysis  systems. 

(2)  Develop  capability  to  support  current  and  future  data  requirements  of 
the  DARPA  research  community. 

(1)  Provide  a  test  bed  for  research  on  more  effective  seismic  analysis 
methods. 

(4)  Meet  the  international  and  special  trilateral  data-exchange  obligations 
anticipated  under  a  Comprehensive  Test  Ban  Treaty  (currently  under 
negotiation). 

(5)  Demonstrate  new  technologies  applicable  to  other  seismic  missions. 

Lincoln  Laboratory  has  been  assigned  the  role  of  overall  technical  leadership  in  the  develop¬ 
ment  of  this  SDC,  which  will  be  carried  out  by  Lincoln  Laboratory  and  other  DARPA  contractors. 

The  accepted  SDC  design  consists  of  a  set  of  minicomputers  interconnected  by  a  local  com¬ 
puter  network.  For  convenience  in  subsequent  discussions,  the  system  can  be  broken  down  into 
seven  subsystems  (Fig.  1-1),  namely: 

(1)  Communications  interface  and  data-logging  subsystem,  which  will  carry 
out  all  input  data  reception,  reformatting,  and  logging  functions. 

(2)  Real-time  analysis  subsystem,  which  will  provide  a  facility  for  automatic 
and  interactive  analysis  of  the  input  data  streams  in  real  time,  if  such  a 
function  is  required. 

(3)  Data  Base  Management  Subsystem,  which  will  provide  both  on-  and 
off-line  storage  of  all  raw  and  processed  data  within  the  system, 
together  with  appropriate  retrieval  mechanisms. 

(4)  Seismic  Analysis  Station,  which  will  be  a  single-user  interactive  wave¬ 
form  and  alphanumeric-data -analysis  system  for  both  waveform  and 
event  processing  functions.  The  final  SDC  will  have  several  of  these 
stations. 

(5)  Local  Computer  Network  subsystem,  which  is  a  high-data-rate  com¬ 
munications  network  linking  all  the  computers  in  the  system. 
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(t)  General  computational  support  subsystem,  which  will  provide  a  variety 
of  computational  services  for  software  development,  certain  tvpes  of 
user  services,  and  svstem  monitoring  and  control. 


(7)  Hesearch  support  subsystem,  which  will  carry  out  all  functions  related 
to  interactive  use  oi  the  svstem  (both  remote  and  in-house)  by  the 
seismologiral  community. 


These  subsystems  are  logical  elements,  defined  in  terms  of  system  functions,  and  should  not  be 
regarded  as  necessarily  representing  specific  hardware  components. 

Lincoln  Laboratory  has  embarked  specifically  on  the  development  of  the  Seismic  Analysis 
Station  (SAS).  the  Local  Computer  Network,  and  the  Data  Base  Management  Subsystem.  Prog¬ 
ress  in  each  of  these  areas  is  described  in  the  sections  that  follow. 

Hardware  for  the  SAS  has  beenseleeted  and  ordered.  and  delivery  is  expected  within  3  months. 
In  the  meantime,  we  have  devoted  our  main  efforts  to  the  development  of  an  interim  system  that 
will  provide  an  event-processing  capability  on  alphanumeric  arrival  data.  We  refer  to  this  in¬ 
terim  system  as  Mod  I  in  the  sections  that  follow'.  It  is  designed,  of  necessity,  to  run  on  our 
in-house  PDP-ll/45  computer,  and  will  provide  a  test  bed  for  the  development  of  data-processing 
techniques  that  will  be  implemented  in  a  more  advanced  form  on  the  prototype  SAS. 

The  Local  Computer  Network  subsystem  has  been  specified  for  procurement  purposes,  and 
industrial  proposals  for  both  hardware  and  software  elements  of  this  subsystem  are  scheduled 
to  be  received  in  April  1980. 

The  data  management  system  design  has  been  refined  to  eliminate  some  of  the  duplication 
of  indexing  information  among  the  subsystems  defined  in  the  special  DARPA  report.  The  effect 
of  this  design  refinement  is  to  simplify  the  interfaces  between  the  subsystems  and  to  simplify 
the  implementation  requirements.  Another  issue  which  has  been  studied  in  this  refinement  is 
efficiency  of  operation  of  the  waveform  data  storage.  This  has  led  to  the  modification  of  the 
proposed  indexing  structure  of  the  waveform  database.  Elements  of  the  parameter  database 
component  of  this  subsystem  are  being  studied  during  the  Mod  I  development  described  above. 
Procurement  procedures  that  may  be  necessary  for  this  subsystem  will  be  initiated  during  the 
remainder  of  the  current  fiscal  year.  A<  Chinnery 

A.  G.  Gann 


B.  THE  MOD  I  INTERIM  SYSTEM:  OVERVIEW 

The  purpose  of  the  Mod  I  system  is  to  provide  a  test  bed  for  the  development  of  the  data 
structures  and  data-processing  routines  that  will  be  implemented  on  the  SAS  later  in  the  current 
fiscal  year  (upon  completion  of  hardware  procurement).  Specifically,  Mod  I  consists  of  a  pre¬ 
prototype  version  of  the  software  that  will  be  required  for  the  storage,  manipulation,  and  pro¬ 
cessing  of  alphanumeric  data.  It  has  been  developed  on  the  in-house  PDP-ll/50  computer  at 
Lincoln  Laboratory's  Applied  Seismology  Group. 

The  basic  analyst  procedures  addressed  by  the  Mod  I  system  are: 

(1)  Construct  an  accessible  database  of  seismic  arrival  parameters, 

(2)  Subject  this  database  to  an  automatic  association  program,  to  produce 
a  preliminary  list  of  events, 
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Provide  for  analyst  review  and  modification  of  these  events,  including 
forming  new  events  from  unassociated  arrivals,  and 

(4)  Output  a  completed  event  bulletin. 

As  a  developmental  system,  one  of  the  main  requirements  on  the  design  of  Mod  I  is  that  it 
should  he  easy  to  add  new*  processing  algorithms  and  procedures.  To  do  this,  we  have  designed 
lata -access  package  that  allows  us  to  easily  change  the  parameters  making  up  the  data  files, 
an  1  we  have  implemented  the  various  processing  steps  as  separate  programs  which  may  be  cas- 
1 1 \  modified  or  replaced.  We  have  made  provision  for  both  automatic  programs  that  run  unat¬ 
tended  on  the  full  database,  and  for  interactive  programs  that  provide  the  operator  with  visual 
feedback  after  each  step.  The  first  version  of  Mod  I  is  minimal  in  the  functions  it  provides, 
but  has  set  the  stage  for  the  easy  addition  of  more  extensive  algorithms  and  procedures. 

The  basic  elements  of  the  Mod  I  database  are  an  arrival  list,  an  event  list,  and  an  associa¬ 
tions  list.  Mach  record  in  the  arrival  list  describes  the  reception  of  one  phase  of  a  seismic 
signal  at  a  single  station.  It  contains  fields  for  the  station,  reported  phase,  arrival  time,  am¬ 
plitude.  period,  etc.  Each  record  in  the  event  list  describes  a  single  event:  its  latitude,  longi¬ 
tude,  origin  time,  depth,  error  estimates,  etc.  Each  record  in  the  associations  list  relates 
one  arrival  to  one  event.  This  includes  association-related  data  such  as  distance,  azimuth, 
magnitude,  and  residual.  A  flag  can  mark  a  given  association  as  defining  (that  is,  used  in 
computing  tne  event  location)  or  nondefining.  The  associations  have  been  placed  in  a  separate 
list,  rather  than  being  incorporated  in  the  arrivals  list  as  is  commonly  done,  to  allow  the  pos¬ 
sibility  of  associating  a  single  arrival  wiui  more  than  one  event. 

Figure  1-2  shews  the  main  software  blocks  that  make  up  the  current  Mod  I  system  and  how 
the  data  flow  between  them. 

Reformat  is  ;.ae  program  that  takes  the  raw  arrival  data  and  transforms  them  into  a  Mod  I 
arrival  list.  The  current  version  of  Reformat  will  handle  two  types  of  test  data  that  we  have 
been  using  to  test  the  Mod  I  system:  data  from  the  International  Seismic  Month  (ISM),  and  syn¬ 
thetic  data  produced  by  a  program  described  in  Sec.  II-C.  After  we  start  receiving  ISDE  data  over 
the  WMO  network.  Reformat  will  be  expanded  to  handle  that  also. 

AA  is  the  Automatic  Association  program.  It  goes  through  the  full  arrival  list  and  attempts 
to  associate  arrivals  to  produce  a  preliminary  event  list  and  associations  list.  The  method  used 
is  described  in  more  detail  in  Sec.  It- A. 

Select  is  a  program  tha.,  given  an  event  number,  searches  the  main  database  for  that  event 
and  all  arrivals  and  associations  relating  to  that  event.  It  outputs  these  as  three  "per-event" 
lists  that  serve  as  the  analyst’s  working  database.  All  further  processing  is  done  on  the  "per- 
event"  lists  until  the  Accept  progrcun  is  called. 

Select  has  several  options.  In  addition  to  selecting  the  arrivals  that  are  currently  associated 
with  the  requested  event,  it  can  also  search  for  and  include  any  unassociated  arrivals  that  have 
an  absolute  travel-time  residual  with  rerpect  to  the  event  of  less  than  a  specified  amount.  It 
can  also  select  only  the  unassociated  arrivals  for  r.  narcicular  time  period,  so  that  the  analyst 
can  then  attempt  to  generate  an  entirely  new  event. 

Display  is  a  program  to  display  and  allow  editing  of  the  "per-event"  lists.  Using  Display, 
the  analyst  can  associate  or  unassociate  any  arrival  with  the  event  under  consideration. 

Locate  is  a  program  to  relocate  the  event,  based  on  any  modifications  to  the  arrivals  or 
associations  that  might  have  been  made  while  in  the  Display  program.  It  has  three  options  for 
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thr  initial-event  location  it  can  use  the  location  currently  in  t^e  database,  it  can  accept  a 
typed -in  initial  location,  or  it  can  use  the  time  and  location  of  ti  first  reported  arrival. 

Accept  is  the  program  which  will  return  the  modified  "per-ever. lists  back  to  the  main 
database. 

Par,  pev,  and  px,  which  are  not  shown  in  Fig.  1-2,  are  programs  to  print  out  the  arrival, 
event,  and  associations  lists. 

Bulletin,  which  is  also  missing  from  the  diagram,  is  a  program  to  print  out  the  database 
in  standard  seismic  bulletin  format.  It  lists  each  event,  followed  by  the  arrivals  associated 
with  that  event.  L.  j.  Turek 


C.  THE  MOD  I  DATA-ACCESS  PACKAGE 

The  programs  in  the  Mod  I  system  use  a  data-access  package  that  has  been  developed  spe¬ 
cifically  for  Mod  I.  These  access  subroutines  are  used  to  read,  write,  ai  t’  search  for  data 
stored  on  disk  in  a  specific  list  format.  Files  in  this  format  are  called  pa>  imeter  files.  In 
the  Mod  I  context,  each  of  the  basic  lists  -  arrival,  event,  and  associations  -  discussed  in  the 
previous  section  is  in  the  form  of  such  a  parameter  file. 

A  feature  of  the  system  is  that  programs  which  work  on  a  set  of  parame  ters  will  work  on 
any  files  that  contain  those  parameters,  regardless  of  what  other  parameters  are  present.  Thus, 
new  oarameters  can  be  added  to  the  database  without  the  necessity  of  recompiling  or  relinking 
existing  programs.  This  flexibility  can  save  much  effort  when  a  system  is  in  its  developmental 
phase. 

When  the  access  subroutine  package  is  used  to  open  a  parameter  file,  it  will  be  given  by  the 
user  a  description  of  the  parameters  of  interest  to  the  program.  It  will  then  set  up  a  mapping 
from  the  parameters  found  in  the  file  to  the  parameters  in  a  structure  provided  by  the  user. 

Once  the  mapping  is  set  up,  all  subsequent  reads  and  writes  to  the  file  are  accomplished  via 
the  user's  structure. 

A  second  level  of  the  access  subroutine  package  is  the  indexing  level,  which  is  constructed 
on  top  of  the  basic  access  level  described  above.  At  the  indexing  level,  the  user  program  may 
choose  to  build  indexes  on  selected  parameters  in  selected  files.  Once  an  index  is  constructed, 
it  may  be  used  to  "look  up"  the  records  in  the  file  that  contain  certain  values.  Functions  pro¬ 
vided  at  this  level  include: 

pirec(EXACT),  which  searches  for  a  record  containing  an  exact  value 
of  the  indexed  parameter. 

pirec(NEAREST),  which  searches  for  the  record  containing  the  nearest 
value  equal  to  or  greater  than  the  given  value  of  the  indexed  parameter. 

pirecs,  which  searches  for  all  records  containing  the  exact  value  of  the 
indexed  parameter. 

pirange,  which  searches  for  all  records  containing  a  value  of  the  indexed 
parameter  within  the  given  range. 

pistep,  which  will  treat  the  file  as  if  ordered  on  the  indexed  parameter 
and  step  forward  or  backward  "n"  records. 
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The  design  of  the  indexing  level  calls  for  implementing  the  indexes  as  B+  trees,11  which  are 
id  cal  Is  suited  for  both  fast  searches  and  the  stepping  function  (which  is  essential  for  ordered 
li.splavs).  Hut  in  order  to  begin  work  on  Mod  I  without  waiting  for  the  indexing  level  to  be  im¬ 
plemented  in  its  final  form,  we  have  provided  for  an  interim  version.  Currently  in  use,  this 
interim  version  has  all  the  functions  that  will  be  available  in  the  final  version,  but  has  imple¬ 
mented  them  without  indexes,  as  brute-force  searches  of  the  parameter  file.  When  the  full 
implementation  is  eventually  achieved,  it  can  be  added  to  the  system  with  no  changes  required 
in  the  calling  programs. 

Although  these  access  routines  were  originally  designed  for  use  with  Mod  I,  they  are  suf¬ 
ficiently  general  that  we  expect  them  to  prove  themselves  useful  in  several  other  contexts.  For 
example,  the  Distributed  Sensor  Networks  project  is  using  them  to  contain  their  waveform, 
spectral,  and  azimuthal  response  data  output  from  their  signal-processing  routines.  We  also 
feel  that  they  could  form  the  basis  for  a  general  data -handling  package  similar  to  the  DADS 
package  that  was  implemented  on  the  PDP-7'sJ 

L.  J.  Turek 
P.  T.  Crames 


D.  LOCAL  COMPUTER  NETWORK 

The  previously  reported  design  for  the  SDC  features  an  architecture  utilizing  multiple  mini¬ 
computers  interconnected  by  a  local  computer  network.  The  local  network  is  needed  to  provide 
an  effective  and  efficient  means  of  transmitting  data  and  control  information  between  the  various 
computers  which  will  make  up  the  SDC.  This  computer  network  must  serve  to  interconnect  the 
computers  very  reliably  and  with  a  high  intercomputer  data  rate.  During  this  reporting  period, 
a  specification  for  procurement  has  been  drawn  up  and  issued.  The  proposals  of  prospective 
vendors  are  due  in  April  1980. 

The  specific  requirements  for  the  multiple-access,  broadcast-mode,  local-computer-packet 
network  needed  for  the  SDC  are: 

(1)  The  network  shall  provide  a  signaling  speed  of  2  Mbps  or  greater.  The 
design  used  shall  be  capable  of  being  extended  to  5  Mbps  or  greater, 
with  no  change  in  host  computer  software. 

(2)  The  design  of  the  interface  shall  be  modular  such  that  interfacing  to  a 
new  host  requires  replacing  only  the  host  specific  part  of  the  interface. 

(3)  The  initial  network  shall  provide  interfaces  for  the  DEC  PDP-11  Unibus 
with  DMA  transfer  of  data,  in  packets  with  maximum  size  of  not  less 
than  512  bytes,  to  and  from  host  memory.  The  initial  installation  shall 
provide  hardware  for  6  host  (Unibus)  connections. 

(4)  The  network  shall  be  capable  of  supporting  30  or  more  hosts.  This  re¬ 
quires  the  ability  to  address  at  least  30  distinct  addresses  on  the  network. 

This  also  requires  that  normal  operation  be  possible  with  a  minimum  of 
30  hosts  physically  attached  to  Uie  communication  medium. 


* D.  Comer,  "The  Ubiquitous  B-Tree/'  ACM  Computing  Surveys  LI,  2  (1979). 

t  Seismic  Discrimination  Semiannual  Technical  Summary,  Lincoln  Laboratory,  M.I.T.  (31  De¬ 
cember  1972),  DDC  AD-757560. 
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(3)  The  network  sha;  tie  to  support  an  "eavesdropping"  mode  of  opera¬ 

tion.  This  means  th.  a  host  can  passively  copy  a  message  addressed 
to  a  specific  address  and  not  perform  the  acknowledgment  protocol. 

This  allows  the  simultaneous  reception  of  data  by  multiple  hosts,  with 
one  host  serving  as  the  primary  recipient  for  protocol  purposes. 

(6)  Additionally,  the  network  shall  support  a  broadcast  mode  where  a 
message  is  directed  to  all  attached  hosts.  The  broadcast  mode  does 
not  require  positive  acknowledgment  of  each  transmission.  Normally 
addressed  transmissions  shall  be  positively  acknowledged  by  the 
addressed  host. 

(7)  The  network  shall  be  capable  of  supposing  physical  separation  between 
the  most  distant  nodes  of  at  least  1  km.  The  initial  installation  will  not 
require  separation  between  the  most-distant  host  to  exceed  30  m. 

(8)  The  interfaces  shall  be  fail-safe,  i.e.,  no  interface  failure  can  disable 
more  than  its  attached  host. 

(9)  The  interface  shall  employ  traffic  sense  access  to  the  communications 
media.  That  is,  the  interface  shall  check  for  the  presence  of  traffic 
before  initiating  any  transmission.  It  should  employ  collision  detection, 
which  is  the  continued  checking  for  interfering  traffic  during  the  period 
of  transmission.  If  it  doesn't  use  collision  detection,  it  should  use  a 
higher  signaling  speed  sufficient  to  duplicate  the  performance  of  the 
design  signaling  speed  (2  Mbps)  using  collision  detection.  If  there  is 
any  doubt  as  to  the  adequacy  of  the  speed  of  a  noncollision  detection 
system,  that  speed  should  be  supported  by  analysis  using  published 
analytic  results  comparing  carrier  sense  and  collision  detection  per¬ 
formance  in  capacity  at  a  fixed  expected  service  delay. 

(10)  Although  it  is  planned  to  utilize  the  network  within  the  confines  of  a 
single  building,  the  interface  shall  employ  isolation  techniques  suf¬ 
ficient  to  allow  attachment  of  hosts  connected  to  completely  separate 
main  power  sources,  i.e.,  different  utility  services  -  ever,  completely 
different  utility  companies. 

(11)  The  physical  design  shall  be  such  that  the  host  interface  can  be  replaced 
within  5  min.  or  switched  to  a  completely  separate  physical  communi¬ 
cation  medium  (previously  installed)  within  5  min.  to  service  a  failure 
in  interface  or  a  communication  medium. 

(12)  The  interface  shall  employ  an  error  detection  and  control  technique  at 
least  as  good  as  the  CRC  technique  specified  in  the  International  Stan¬ 
dards  Organization's  IIDLC  (ISO  document  2309). 

(13)  The  interface  shall  have  provision  for  fair  access  to  the  local  computer 
net  services  in  the  event  of  extended  saturation  of  load. 

(14)  The  interface  should  allow  the  incorporation  of  a  priority  scheme  to  allow 
favored  access  to  the  communication  medium  by  some  host  in  overload 
conditions. 
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Tlu*  prefer  rod  design  of  the  local  computer  network  to  satisfy  these  requirements  is  a  set 
of  interlaces  to  a  passive  communication  medium.  Each  interface,  one  lor  each  host  computer 
on  the  net,  would  implement  the  low-level  protocols  to  transfer  messages  between  hosts  with  a 
minimum  of  host  processing.  The  host  software  would  implement  higher  level  protocols  to  ex¬ 
change  requests  and  responses  between  hosts  and  to  transfer  logical  collections  of  data  between 
hosts  using  the  lower-level,  hardware-implemented  protocols. 

Unfortunately,  there  are  no  commercially  available  local  computer  networks  that  meet  the 
desired  characteristics.  The  procurement  currently  under  way  is  intended  to  identify  the  best 
compromise  of  cost  and  performance,  taking  into  account  the  required  host  software  as  wrell  as 
the  hardware  costs,  that  will  satisfy  these  requii  ements.  As  soon  as  the  hardware  is  chosen, 
the  design  of  the  protocols  and  the  host  software  to  completely  realize  the  needed  intercomputer 
communication  capability  will  be  started.  A.  G.  Gann 

E.  DESIGN  FOR  THE  INITIAL  DATABASE  PROTOTYPE 

Database  design  efforts  have  focused  on  the  initial  prototype  database  implementation.  Pre¬ 
sented  below  is  an  outline  of  the  design  and  a  discussion  of  the  role  of  this  prototype  in  the  over¬ 
all  plan  of  database  development.  The  initial  database  module  will  be  implemented  as  a  series 
of  functional  units  whose  interface  will  allow  each  individual  unit  or  a  selected  set  of  units  to  be 
isolated  on  its  own  computer.  The  system  will  be  implemented  on  a  modified  UNIX.  Initially 
it  will  be  co-resident  with  the  SAS(  but  should  be  capable  of  being  separated  therefrom  at  an 
early  date. 

The  goals  of  the  initial  prototype  implementation  are  to  provide  a  useful,  nonredundant  im¬ 
plementation  for  use  by  other  prototype  elements  of  the  SDC  ana  to  get  information  which  will  be 
useful  in  building  high-performance  fully  redundant  implementations  which  still  provide  the 
same  basic  functions.  Ir  order  to  fulfill  this  second  goal,  the  initial  implementation  will  attempt 
as  fine-grained  a  division  into  functional  units  as  possible.  This  will  allow  freedom  in  packaging 
machine  functions  in  later  implementations.  Isolated  functional  units  will  allow  multiple  active 
copies  of  a  function  to  reside  in  the  same  host.  This  will  allow  numerous  experiments  to  be 
conducted  to  determine  the  suitability  of  many  proposed  designs  for  fully  redundant  operation, 
before  implementing  the  final  version. 

1.  Functional  Design 

Functional  units  are  isolated  in  such  a  fashion  that  they  may  be  placed  on  separate  hosts  when 
connected  by  a  local  networr.  Functional  units  communicate  by  means  of  messages  which  are 
passed  between  units  by  the  kernel,  in  the  same  fashion  that  the  local  network  passes  messages 
between  hosts. 

Co-residence  of  multiple  functional  units  on  a  single  host  will  provide  a  fruitful  environment 
for  testing  the  system's  response  to  unit  or  communications  failure.  Provision  of  commands  to 
terminate  or  restart  a  given  unit  or  disable  the  communication  medium  will  allow  the  robustness 
of  the  system  to  be  established  for  large  configurations  which  could  not  easily  be  realized  if  a 
host  were  necessary  for  each  functional  unit. 

Performance  predictions  can  also  be  supplied  by  this  type  of  environment  by  using  virtual- 
time-management  techniques  similar  to  those  used  by  the  virtual  machine  emulator.’  This  will 
allow  determination  of  performance  under  a  wide  range  of  assumptions. 

*  M.  D.  Cannon  et  aL,  "A  Virtual  Machine  Emulator  for  Performance  Evaluation,"  Commun. 

ACM  23,  71  -80  (1980). 
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2.  Waveform  Database 


The  Waveform  Database  subsystem  is  one  of  two  sets  of  intercommunicating  functional  units, 
the  other  being  tne  Parameter  Database  subsystem.  The  Waveform  Database  is  responsible  for 
the  storage  of  all  wa\  jform  data  within  the  system,  including  reference  event  data.  Ail  users 
of  the  Waveform  Database  will  direct  their  requests  to  a  single  point  of  interface,  the  index 
manager.  The  index  manager  is  responsible  for  selecting  a  disk  or  tape  manager  which  has 
the  data  and  which  will  satisfy  the  request. 

3.  Index  Manager 

The  index  manager  maintains  an  index  of  all  waveform  data  in  the  system,  whether  stored 
on  disk  or  tape.  This  allows  it  to  designate  the  appropriate  unit  to  respond  to  a  waveform  read 
request.  The  selected  unit  can  then  initiate  a  connection  to  the  requestor  and  transmit  the  data 
independently. 

The  index  manager  is  also  responsible  for  the  allocation  of  disk  space.  The  index  manager's 
control  of  disk -space  allocation  makes  it  responsible  for  fielding  waveform  write  requests  from 
the  communications  interface.  Upon  receiving  a  write  request,  the  index  manager  allocates  a 
disk  extent,  thereby  selecting  the  disk  manager  responsible  for  the  spindle  on  which  the  extent 
is  to  reside.  A  connection  is  then  set  up  and  data  transfer  proceeds  thereon. 

The  index  manager's  control  of  disk  allocation  extends  management  of  the  bin  tape  crea¬ 
tion  process.  Data  write  requests  for  each  data  extent  pass  through  the  index  manager,  and  it 
is  therefore  in  a  position  to  maintain  a  queue  of  data  to  be  dumped  and  to  present  a  manifest 
containing  the  location  of  such  data  to  the  tape  module  so  that  it  may  produce  a  bin  tape.  Once 
data  have  been  so  dumped,  the  index  manager  will  keep  track  of  the  use  of  each  extent  and  may 
thereby  deallocate  extents  on  a  least  recently  used  (LRU)  basis. 

4.  Disk  Manager 

The  disk  manager  is  responsible  for  the  high-speed  transfer  of  waveforms  upon  request  by 
the  index  manager.  The  disk  manager  sets  up  a  connection  with  the  original  requestor  identified 
by  the  index  manager  and  then  proceeds  to  transfer  the  data  between  the  established  network 
connection  and  the  disk.  Data  transfer  proceeds  sequentially  within  a  single  extent,  so  no  in¬ 
dexing  is  needed.  At  the  completion  of  the  transfer,  the  index  manager  is  notified  so  that  fur¬ 
ther  segments  may  be  transferred  in  a  like  manner. 

5.  Tape  Manager 

The  tape  manager  is  responsible  for  producing  and  updating  bin  tapes  in  response  to  index 
manager  instructions,  and  for  reading  data  on  bin  tapes  and  reference  event  tapes  back  into  the 
system.  It  is  assumed  here  that  other  use  of  tapes  will  be  handled  in  a  logically  distinct  func¬ 
tional  unit.  This  may  be  on  another  machine  or  possibly  as  a  separate  task  on  the  same  machine. 

6.  Parameter  Database 

The  Parameter  Database  subsystem  is  responsiole  for  system  storage  of  nonwaveform  data. 
It  is  layered  into  two  functional  units:  the  file  system,  which  provides  a  UNEX-type  file  system 
as  a  central  storage  repository  for  network  users;  and  the  database  manager,  which  uses  the 
underlying  file  system  and  provides  an  object-oriented  data-management  system. 
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The  file  system  will  be  based  on  the  current  UNIX  kernel  file  system  code.  By  moving  this 
{inu  tinu  into  a  supervisor  address  space  and  making  it  a  server  responding  to  messages  from 
the  other  units  passed  to  it  by  the  kernel,  the  basis  is  laid  for  a  central  file  repository  for  the 
local  network.  This  repository  may  be  used  to  read  or  write  individual  file  records  or  to  copy 
entire?  files  from  their  canonical  version  to  an  SAS  or  other  network  user. 

The  database  manager  provides  an  object-oriented  data -independent  interface,  to  be  used 
primarily  in  the  handling  of  event  and  arrival  lists.  The  searching  of  the  database  for  objects 
which  have  certain  characteristics  will  be  performed  by  the  database  manager  upon  command, 
with  the  matching  objects  returned  to  the  caller  in  a  caller-specified  format.  The  precise  re¬ 
quirements  of  the  application  in  this  regard  have  yet  to  be  determined. 

D.  B.  Noveck 


F.  TUB  SDC  DOCUMENTATION  SYSTEM 

During  the  initial  design  phase  of  the  SDC,  it  was  perceived  that  an  organic  procedure  was 
required  to  assure  proper  information  flow  among  the  designers  of  the  center.  This  is  espe¬ 
cially  important  since  some  of  the  designers  and  implementors  of  sub-projects  are  not  located 
at  the  Applied  Seismology  Group.  Also,  a  method  of  exercising  a  relatively  tight  overview  of 
the  project's  progress  was  desired.  An  accessible  documentation  system  was  seen  as  an  im¬ 
portant  factor  in  achieving  this  objective. 

The  aim  of  the  documentation  system  is  to  have  a  flexible  system  which  provides  a  "friendly 
face"  to  both  the  generator  and  reader  of  the  documents.  It  should  provide  the  means  for  easy 
transfer  of  information  between  all  personnel  in  the  project  throughout  the  history  of  the  project. 
We  decided  to  maintain  the  design,  scheduling,  progress,  and  operations  documentation  on-line, 
on  the  Group's  general-purpose  computers.  Nonresident  personnel  would  then  be  able  to  access 
and  create  documents  via  the  ARPAnet.  These  computers,  of  the  PDP-11  family,  use  the  UNIX 
operating  system  which  has  been  utilized  in  the  past  for  the  implementation  of  text -handling 
systems  (for  example,  AT&T  has  used  the  system  for  publication  of  the  Bell  System  Technical 
Journal). 

The  design  of  the  documentation  system  falls  into  three  main  categories: 

(1)  A  storage  mechanism  for  the  documents  which  is  organic  to  the  design 
of  the  SDC.  The  UNIX  file  system  was  easily  tailored  to  perform  this 
task.  A  documentation  hierarchical  file  structure  was  constructed  which 
parallels  the  overall  system  design. 

(2)  A  set  of  commands  to  easily  access  and  maintain  the  documentation.  A 
level  of  protection  was  required  so  that  documents  which  numerous 
people  might  have  to  modify  would  not  become  cloned,  and  then  modified 
in  parallel;  this  would  result  in  multiple  "up-to-date"  versions  of  the 
same  document. 

(3)  A  method  of  implementing  a  set  of  standards  so  that  the  documents  of 
differing  authors  would  have  a  uniform  appearance.  It  was  decided  to 
use  the  nroff  program,  a  relatively  sophisticated  text  assembly  program, 
to  aid  in  achieving  this  task. 
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1.  Tin*  Hierarchy  of  tin*  Kile  Structure 


The  file  structure  parallels  the  functional  structure  of  the  proposed  SDC.  The  top-level 
directory  is  divided  into  eight  sub-directories,  six  of  which  are  related  to  the  functional  sepa¬ 
ration  of  the  SDC'.  lOach  of  these  contains  sub-directories,  and  not  documents.  Of  the  two  re¬ 
maining,  one  is  for  documents  relating  to  the  SDC,  as  a  whole,  and  the  other  is  a  "catch  all" 
category  for  random  notes  which,  when  generated,  didn't  seem  to  fit  anywhere.  The  eight  di¬ 
rectories  are: 

oview  executive  summary,  system  status,  milestones 

comm  communications  interface 

vsrf  VSC  seismic  research  facility 

sas  seismic  analysis  station 

dbs  database  system 

net  local  network 

gsc  Lincoln  Laboratory  seismic  research  facility 

notes  system -w'ide  information  not  in  any  other  category 

Within  each  of  the  six  functional  task  directories,  there  are  five  sub-directories: 

oview 

if 

fen 

dgn 

notes 

For  each  task,  files  in  the  "oview"  directory  store  the  current  status  and  milestone  scheduh 
for  that  task  as  well  as  a  brief  summary  of  the  purpose  and  design. 

The  "if"  directory  contains  one  or  more  documents  which  (together)  constitute  the  interface 
specifications  for  the  interfaces  to  other  elements  in  the  system.  The  "fen"  directory  contains 
files  which  constitute  the  functional  specification  of  the  overall  subsystem  and  each  separately 
identified  subdivision  of  the  subsystem. 

The  "dgn"  directory  contains  files  documenting  the  design  of  each  part  of  the  subsystem. 
These  documents  give  the  higher  level  design  considerations  for  the  programs.  The  "notes" 
di rectory  is  used  (at  all  levels)  to  catch  thoughts  and  ideas  which  are  not  readily  categorized 
into  the  above  cat* ,  ories.  Some  of  these  may  be  moved  later  into  the  regular  categories,  others 
may  cause  new  ce.egories  to  be  added,  and  others  may  be  dropped  after  discussion  and  further 
thought.  The  purpose  of  this  "catch-all"  category  is  to  ensure  that  no  idea  is  lost  for  lack  of 
a  category. 

It  is  anticipated  that  each  directory  will  have  a  person  who  is  responsible  for  maintaining 
the  integrity  of  the  directory.  This  is  primarily  ensuring  that  all  files  in  the  directory  are 
consistent  and  properly  maintained. 

2.  Documentation  Maintenance  Commands 

The  following  is  a  quick  overview  of  the  commands  used  to  access  and  maintain  the  docu¬ 
mentation  in  this  system.  The  commands  are  reserve,  distribute,  runoff  and  look,  and  since. 
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Reserve 


Makes  a  copy  of  the  named  file  available  for  exclusive 
access  by  a  user.  The  previous  version  of  the  file  (if 
any)  remains  accessible  by  all  users  for  viewing  only. 

Installs  the  new  version  of  a  file  for  access  by  all  users. 

The  previous  version  is  archived  in  a  special  directory 
for  later  copy  to  magnetic  tape.  Thus,  a  complete  audit 
trail  of  versions  of  a  file  is  available. 

Reformats  a  file  from  the  internal  format  to  the  display 
format.  The  output  is  via  the  UNIX  Standard  Output, 
hence  is  available  for  processing  by  filters  or  filing  via 
shell  options.  A  special  variant  of  runoff,  called  look, 
is  used  if  and  only  if  the  file  is  currently  reserved  to 
the  user,  while  runoff  will  give  the  publicly  accessible 
version  at  all  times. 

Gives  a  list  of  all  files  changed  in  the  specified  number 
of  days.  A  specified  directory  or  the  current  directory 
and  all  subsidiary  directories  are  scanned  to  determine 
all  files  which  have  been  changed  (or  created)  within 
the  past  specified  number  of  days. 

3.  Document  Production 

In  order  to  have  a  uniform  appearance  for  the  SDC  documentation,  we  will  use  a  standard 
documentation  method.  This  will  enable  us  to  distribute  up-to-date  versions  which,  if  we  main¬ 
tain  the  documentation  database,  will  accurately  and  clearly  reflect  the  current  status  and  future 
of  the  SDC.  In  addition,  it  will  enable  us,  by  simply  joining  relevant  sections,  to  produce  what¬ 
ever  reports  are  required. 

An  expository  flexible  format  is  required  for  transmitting  conceptual  information.  The 
author's  mind  view  of  the  concept  has  much  to  do  with  the  way  in  which  it  should  be  presented. 
Minimizing  the  transformation  from  the  creator's  mind  to  the  page  will  tend  to  ease  the  docu¬ 
mentation  process  for  the  author,  and  allow  the  concepts  to  be  more  readily  grasped.  Using  a 
common  documentation  package,  which  gives  us  the  flexibility  required,  will  both  ease  the  tran¬ 
sition  from  mind  to  media  and  provide  a  uniform  appearance  to  the  documentation. 

We  will  use  the  nroff  program  which  is  controlled,  from  within  the  text  stream,  by  a  detailed 
group  of  commands.  These  commands  provide  the  means  of  formatting  the  document  to  be  dis¬ 
played.  The  commands  can  be  grouped  together  into  interactive  packages  called  macros. 

We  use  a  subset  of  the  standard  -ms  macro  package,  to  which  additions  have  been  made  to 
fit  our  needs.  Using  a  common  documentation  package,  and  a  few  guidelines,  gives  the  documen¬ 
tation  a  standard  appearance  without  being  a  burden  on  the  authors,  and  gives  the  flexibility 
required. 


Distribute 


Runoff 


Since 


4.  Initial  Results 

The  initial  use  of  the  documentation  system  has  been  very  favorable.  It  was  first  used  to 
document  the  documentation  system.  This  "shakedown"  process  resulted  in  modifications  to 
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components  of  the  system  and,  hence,  to  concomitant  revision  of  the  documentation  system's 
documentation.  The  system  is  now  stable.  Existing  SDC  documentation,  which  was  developed 
outside  the  documentation  system,  is  now  being  moved  onto  the  system. 

R.  S.  Blumberg 
A.  G.  Gann 


Fig.I-1.  Functional  subsystems  of  SDC.  These  may  not  correspond  to  hardware 
components,  since  in  some  possible  configurations  several  functions  may  be  carried 
out  on  same  computer. 


Fig.  1-2.  Software  elements  of  Mod  I  interim  system. 
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II.  SEISMIC'  ALGORITHMS 


A.  AUTOMATIC  ASSOCIATION 

One  of  the  main  goals  of  the  SDC  is  to  produce  a  list  of  earthquake  locations  from  a  list  of 
seismic-wave  arrival  times.  This  process  is  known  as  "association,"  since  sets  of  incoming 
arrivals  must  be  associated  with  particular  epicenters.  About  five  different  automatic-association 
methods  have  been  developed  by  various  seismic  institutions.  Unfortunately,  none  of  the  methods 
has  established  their  superiority  over  the  others,  and  automatic  association  is  thus  an  area  of 
active  research.  The  following  is  a  progress  report  on  the  automatic-association  method  cur¬ 
rently  being  investigated  by  the  Applied  Seismology  Group  at  Lincoln  Laboratory. 

One  of  the  major  problems  in  automatic  association  is  due  to  the  fact  that  earthquakes  may 
occur  at  nearly  the  same  time  anywhere  in  the  world.  For  example,  for  the  International  Seis- 
mological  Month  (ISM)  about  30  percent  of  the  events  overlap  in  time.  Thus  besides  locating  the 
earthquakes  from  the  arrival  times,  which  by  itself  is  not  a  simple  task,1  one  must  decide  which 
set  of  arrival*  should  be  assigned  to  which  earthquakes. 

The  association  problem  can  be  represented  as  the  following  combinatorial  optimization 
problem:  Find  X,  the  set  of  earthquake  locations,  which  minimizes 

Z  +«|X|  +0|U| 

x.eX 

where  A.  is  the  set  of  arrivals  associated  with  location  x.,  P(x.,  A.)  is  the  cost  of  associating 
3  3  3  3 

location  x.  with  the  set  of  arrivals  A.,  a  and  |5  are  appropriate  constants,  and  U  is  the  set  of 

J  J 

unassociated  arrivals.  The  second  and  third  terms  essentially  apply  Occam's  razor  by  trying 
to  associate  as  many  arrivals  as  possible,  and  at  the  same  time  keeping  the  number  of  events 
as  small  as  possible.  In  words,  this  equation  says  v/e  want  to  find  the  smallest  number  of  earth¬ 
quakes  which  best  explains  as  many  arrivals  as  possible. 

Combinatorial  optimization  problems  are  usually  quite  difficult  to  solve,  and  this  one  has 
the  additional  complication  of  errorful  data.  Normally,  suboptimal  heuristic  solutions  are  found 
using  a  form  of  the  hypothesize-and-test  paradigm.  A  small  subset  of  arrivals  (3  to  5)  or  one 
or  more  array  measurements  is  used  to  determine  a  hypothetical  epicenter.  If  enough  other 
arrivals  agree  with  this  epicenter,  it  is  declared  an  event  and  the  arrivals  are  associated  with 
it.  This  process  is  then  repeated  for  any  remaining  unassociated  arrivals.  Although  this  pro¬ 
cedure  is  suboptimal,  it  appears  to  produce  about  80  percent  of  the  events  that  an  analyst  would 
produce  with  an  acceptable  false-event  rate  (G.  Dumphy  of  USGS,  personal  communication). 

1.  The  Method 

The  initial  version  of  the  automatic-association  program  (AA)  was  designed  to  utilize  only 
P  and  PKP  arrival  times.  This  simplifies  the  design  somewhat,  but  still  provides  a  useful  pro¬ 
gram  (90  percent  of  the  ISM  arrivals  are  P  or  PKP  arrivals).  Additional  information  such  as 
amplitude,  later  phases,  and  array  measurements  can  be  incorporated  into  the  algorithm  later 
on  without  much  difficulty. 

The  AA  takes  two  lists  as  input  -  a  list  of  arrivals  and,  optionally,  a  preliminary  list  of 
events.  The  output  of  AA  is  a  list  of  events  and  their  associated  arrivals. 
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If  a  preliminary  event  li.  t  is  available,  arrivals  are  matched  to  events  in  the  list  before  new 
event*  art'  searched  for.  In  this  way,  any  external  event  information  can  be  used,  or  any  arrivals 
missed  during  a  previous  AA  run  can  be  assigned  to  the  appropriate  event  during  a  later  AA  run. 

When  searching  for  an  event,  the  AA  inspects  each  unassociated  arrival  in  turn.  If  this 
arrival  (called  the  key  arrival)  is  consistent  with  enough  other  arrivals,  it  is  used  to  guide  a 
search  for  an  event  (see  below).  Arrivals  are  temporarily  allowed  to  be  associated  with  any 
number  of  events.  Later  on,  an  "event  critic"  stage  is  used  to  remove  multiple  associations 
and  redundant  events. 

The  event  critic  inspects  each  arrival  and  assigns  it  only  to  the  largest  event  with  which  it 
associates.  Any  events  which  no  longer  have  enough  associated  arrivals  are  discarded.  Thus, 
the  cx’itic  acts  as  a  form  of  Occam's  razor.  Typically,  a  discarded  event  is  caused  by  a  false 
alarm  which  fortuitously  combines  with  several  real  arrivals  to  produce  a  bogus  event.  This 
critic  step  is  vital  because  it  provides  A  A  with  a  more  "global"  viewpoint  than  it  would  get  from 
a  simple  left-to-right  pass  over  the  arrival  list.  Without  it,  many  false  events  can  be  produced 
and  true  events  missed. 

To  search  for  an  event,  a  circular  grid  pattern  of  epicenters  is  placed  around  the  location 
of  the  station  corresponding  to  the  key  arrival  as  shown  in  Fig.  II- 1.  The  origin  time  at  each 
epicenter  is  chosen  as  if  the  key  arrival  had  originated  from  that  spot.  The  arrivals  are  then 
matched  against  each  hypothetical  epicenter.  An  arrival  is  associated  with  an  epicent<  •  if  it 
matches  it  within  a  specified  tolerance  which  depends  on  the  grid  spacing.  The  epicenter  which 
agrees  with  the  most  arrivals  (or  has  the  lowest  cost  in  case  of  a  tie)  is  used  as  the  center  of 
further  ciicular  search  with  a  reduced  grid  spacing.  The  process  is  repeated  until  either  too 
few  arrivals  remain  associated  or  a  sufficiently  accurate  epicenter  is  obtained.  This  technique 
is  similar  to  one  used  by  Snell. 

One  advantage  of  this  technique  is  that  it  simultaneously  associates  arrivals  to  an  event  and 
produces  a  location  for  the  event.  Also,  it  avoids  some  of  the  problems  found  by  traditional  loca¬ 
tion  methods  when  using  errorful  arrival  data. 

2.  Results 

The  testing  of  association  algorithms  must  be  done  with  care  since  the  test  data  can  exert 
a  strong  and  unnoticed  influence  on  the  development  of  the  algorithm.  The  AA  will  be  tested  on 
both  actual  seismic  data  from  the  ISM  and  synthetic  seismic  data.  The  ISM  data  set  is  partic¬ 
ularly  valuable  because  it  has  been  carefully  associated  by  seismic  analysts  and  contains  over 
40,000  arrivals  from  nearly  1,000  events. 

The  preliminary  results  on  a  little  over  800  ISM  arrivals  is  quite  encouraging.  Table  II- 1 
compares  how  arrivals  are  interpreted  by  the  ISM  analysts  and  by  AA.  For  example,  the  table 


TABLE  ll-l 

APPLICATION  TO  ISM  ARRIVAL  DATA. 
COMPARISON  OF  AUTOMATIC  AND  MANUAL 
ARRIVAL  ASSOCIATIONS  (PERCENT) 


AA 

Unassociated 

Associated 

Unassociated 

30.1 

9.7 

Associated 

5.1 

55.2 
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TABLE  11-2 


COMPARISON  OF  EVENTS  FOUND  BY  AA 
WITH  THOSE  FOUND  BY  THE  ISM  ANALYSTS 


Number 

Percent 

MATCH 

20 

62.5 

NEW 

5 

15.6 

SPLIT 

7 

21.9 

MERGE 

4 

12.5 

MISS 

1 

3.1 

shows  that  5.1  percent  of  the  arrivals  associated  with  events  by  ISM  analysts  remained  unasso¬ 
ciated  after  processing  by  AA. 

Table  II- 2  compares  the  events  found  by  ISM  analysts  and  by  AA.  The  terms  in  the  left- 
hand  column  have  the  following  interpretation: 

MATCH  The  set  of  associated  arrivals  were  either  unassociated 

by  ISM  analysts  or  were  associated  to  the  same  event  by 
ISM  analysts. 

NEW  The  associated  arrivals  were  mostly  unassociated  by  ISM 

analysts. 

SPLIT  Ac  rivals  from  one  ISM  event  were  split  into  multiple  events 

by  AA. 

MERGE  The  et  *<*  associated  arrivals  came  from  more  than  one 

ISM  <. vent. 

MISS  The  ISM  event  was  completely  missed  by  AA. 

These  preliminary  results  are  encouraging;  only  one  small  event  (six  arrivals)  was  missed 
by  AA  and  over  five  additional  events  were  found.  Three  of  these  events  had  eight  or  more  as- 

sociated  arrivals.  K.  R.  Anderson 

N.  S.  Snell 

B.  AN  IMPROVED  DETECTOR  FOR  EMERGENT  ARRIVALS 

Most  proposed  detectors  employ  the  Z  statistic  to  compare  the  average  power  over  a  short 
window  (~1. 5  sec)  against  a  threshold.  To  make  selection  of  this  threshold  easier,  the  short¬ 
term  average  power  (STAP)  is  normalized  to  Z  using  Eq.  (II- 1);3 

Z  =  (STAP  -  MUstap)/SIGMAstap  .  (II- 1) 

Therefore,  Z  has  a  mean  of  zero  and  a  standard  deviation  of  one.  MU  and  SIGMA  are  the 
long-term  properties  of  the  STAP  and  are  given  initial  value  during  detector  startup.  They  are 
updated  using  Eqs.  (II- 2)  and  (II-  3) :3 


MUn+1  =  (1  -  E)  *  MUn  +  E  *  STAPn 
SIGMAn+1  =  (1  —  E)  *  SIGMAn  +  E  *  (STAPn  -  MUn)  **  2 
For  a  given  startup  population  N,  E  =  l/N. 


(H-2) 
(II- 3) 


Our  detector,  as  well  as  the  SRO  detector,  is  a  Z  statistic  detector.  In  order  to  accurately 
test  the  performance  of  our  detector  we  acquired  6  h  of  continuous  broadband  data  from  the  A1 
huquerque,  Nt>w  Mexico  SHU  (ANMO).  These  data  were  filtered  to  exactly  simulate  a  short- 
period  (SP)  record.  The  continuous  SP  record  was  then  3-pole  Butterworth  bandpass  filtered 
from  0.5  to  2.0  Hz,  the  same  band  used  by  the  SRO  detector. 

We  now  had  a  way  to  directly  compare  the  SRO  detector  with  other  detectors  (human  or  ma¬ 
chine).  Human  analysis  of  the  6-h  record  showed  that  the  SRO  detector  was  unable  to  pick  up  all 
emergent  signals.  These  emergent  signals  are  not  small.  The  problem  is  that  slowly  emergent 
signals  appear  to  the  detector  as  would  a  slow  increase  in  noise.  With  each  slightly  higher  value 
for  the  STAP,  the  values  of  MU  and  SIC.MA  are  pulled  up.  Eventually  the  whole  signal  is  passed 
over,  MU  and  SIGMA  are  now  contaminated  with  signal,  and  the  detector  is  less  sensitive  for  a 
while.  To  trap  these  emergent  signals,  we  introduced  a  time  lag  between  the  STAP  being  tested 
and  the  STAP  used  to  update  Ml!  and  SIGMA. 

Figure  II- 2  shows  that,  with  the  SRO  detector  at  point  B,  MU  and  SIGMA  are  already  con¬ 
taminated  with  signal.  This  being  the  case,  there  is  little  chance  of  a  detection  by  point  C.  With 
the  lag  introduced,  when  the  detector  is  at  point  C  the  values  of  MU  and  SIGMA  have  only  been 
updated  to  point  A.  This  allows  detection  of  emergent  signals  as  small  as  the  smallest  impulsive 
signals  detected.  The  time  lag  used  should  not  be  too  long  since  noise  levels  can  rapidly  increase 
over  a  short  period  of  time,  as  with  the  passage  of  atmospheric  fronts.  By  trial  and  error,  we 
found  the  optimum  value  of  the  lag  to  be  13  sec. 

An  example  of  an  emergent  signal  that  the  SRO  detector  missed,  which  was  picked  up  by 
our  detector,  is  given  in  Fig.  II- 3. 


It  is  worth  noting  that  this  increased  sensitivity  to  emergent  signals  was  not  accompanied 


by  an  increase  in  false  alarms. 


M.  A.  Tiberio 


R.G.  North 


C.  GENERATION  OF  A  SYNTHETIC  ARRIVALS  LIST 

Two  programs  have  been  written  —  the  first  of  which  generates  a  synthetic  list  of  events, 
and  the  second  of  which  computes  the  arrivals  from  these  events  seen  at  a  prescribed  network 
of  stations.  The  purpose  of  these  programs  is  to  provide  a  systematic  method  for  the  evalua¬ 
tion  of  association  and  location  algorithms.  Programs  such  as  this  have  been  written  before, 
but  arc  described  here  since  the  choice  of  parameters  used  is  rather  subjective. 

The  synthetic  events  list  produces  a  set  of  events  with  approximately  correct  geographical 
distribution,  depth  distribution,  and  frequency -magnitude  characteristics.  The  latter  are  com- 

4 

puted  using  the  frequency-moment  relation  proposed  by  Chinnery  and  North,  and  the  moment- 
magnitude  relations  of  Chinnery.5  M  and  m  are  computed  from  the  moment  value  stochasti- 
cally,  to  allow  for  a  range  of  possible  source  mechanisms.  The  number  of  events  listed  during 
a  given  time  period  (say,  1  day)  is  controlled  by  a  threshold  moment  value.  Random  numbers 
are  used  extensively,  so  that  each  run  of  the  program  generates  a  different  event  list. 

One  such  event  list  for  a  single  day  is  illustrated  in  Figs.  II- 4  and  II- 5.  Figure  II- 4  shows 
the  geographical  distribution  of  the  226  events  on  this  list.  The  moment  threshold  chosen  was 
20.0  dyn-cm,  which  corresponds  roughly  to  =  3.3.  The  frequency-m^  plot  for  these  events 
is  shown  in  Fig.  II- 5.  The  "b-value11  of  0.95  is  very  typical  for  real  events.  In  this  list,  there 
were  3  events  with  5,  and  36  events  with  ^  4.  This  corresponds  to  about  13,000  events 

per  year  with  m^  ^  4,  which  is  probably  a  reasonable  estimate.^ 


16 


The  arrivals  program  computes  P-wave  arrivals,  and  some  later  phases,  at  a  network  of 

stations.  For  simulation  purposes,  the  network  used  for  testing  consists  of  50  globally  distrib- 

7  8 
uted  real  stations.  Station  amplitude  biases  (North  )  and  travel- time  anomalies  (Dziewonski  ) 

are  included.  Each  station  is  assigned  an  average  noise  level,  and  (in  the  present  version  of  the 
,  rogram)  is  assumed  to  have  an  infinite  dynamic  range.  The  amplitude  of  each  event  at  each 
station  is  computed,  including  a  Gaussian  component  to  simulate  scattering,  and  an  arrival  is 
declared  when  the  signal  amnlitude  is  substantially  above  the  noise  level.  Travel  times  are  com¬ 
puted  from  standard  tables,  with  a  Gaussian  random  component  added  to  simulate  regional  vai'i 
ations.  A  variety  of  additional  complications  include  removing  each  station  from  service  for 
1  h  per  day,  and  occasional  operator  errors  in  both  depth  phase  identification  and  timing. 
Recause  of  the  random  quantities,  each  arrivals  list  from  a  given  events  list  is  different. 

A  typical  arrivals  list  produced  from  the  226  event  list  mentioned  above  contained  454 
p  arrivals  and  69  later  phases.  Seventeen  events  had  at  least  5  arrivals,  and  could  therefore  (in 
principle)  be  associated  into  events.  Approximately  25  percent  of  the  arrivals  were  from  events 
with  less  than  5  arrivals  and,  in  practice,  would  presumably  be  unassociated.  Smaller  local 
events,  not  included  in  the  event  list,  would  presumably  increase  this  to  about  50-percent  unas- 
sociated  arrivals,  which  is  commonly  observed.  One  immediate  conclusion  is  ihat  identification 
and  removal  of  small  local  events  should,  if  possible,  be  carried  out  before  the  process  of  event 
association. 

Efforts  to  make  the  final  arrivals  list  more  realistic  are  continuing.  Swarms  and  aftershock 
sequences  will  be  added,  and  station  noise  levels  and  dynamic  ranges  refined.  The  programs 
will  be  used  to  test  automatic- association  algorithms,  and  to  assess  the  error  bounds  generated 
by  location  programs.  They  will  also  be  used  to  examine  the  characteristics  of  existing  catalogs 

C 

such  as  that  issued  by  the  International  Seismological  Center.-' 

M.A.  Chinnery 
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Fig.  II- 1.  In  the  association  algorithm  discussed  in  text,  after  a  key 
arrival  is  identified  at  a  particular  station,  a  grid  of  trial  epicenters 
around  that  station  is  used  to  search  for  a  best  event  solution.  Trial 
epicenters  are  arranged  in  a  circular  pattern  as  shown  here.  Typical 
circle  spacings  are  5°  to  10°. 


Fig.  II- 2.  Schematic  of  emergent  arrival. 
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Fig,  II- 3.  Example  of  emergent  arrival  at  AN  MO  which  was  detected  only 
after  introduction  of  lag  in  determination  of  MUstap,  SIGMAstap. 
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III.  SEISMIC  SOURCE  CHARACTERIZATION 


ANOMALOUS  SURFACE-WAVE  RADIATION 
I'ROM  EASTERN  KAZAKH  EXPLOSIONS 

Some  underground  explosions  produce  surface -wave  radiation  which  is  radically  different 
from  *V't  for  others  in  the  same  area.  Such  differences  include  phase  reversal  of  the  Rayleigh 
waves,1 2 *’4'  enhanced  Love-wave  generation,  apparent  time  delays  of  the  Rayleigh  waves  by  a  few 

>  -j  4 

seconds,  '  and  sharp  azimuthal  dependence  of  Rayleigh-wave  amplitudes.  Explanations  for 
these  variations  include  spall  closure  for  the  phase  reversals  and  time  delays,^  and  "tectonic 

A 

strain  release"  for  the  Love-wave  enhancement  and  Rayleigh  azimuthal  variation. 

We  present  here  a  detailed  description  of  differences  observed  between  presumed  nuclear 
explosions  within  the  eastern  section  of  the  Eastern  Kazakh  test  area.  These  events  all  take 
place  within  an  extremely  restricted  region  (improved  locations  for  these  events  are  described 
elsewhere  in  this  report)  about  15  km  square.  We  have  obtained  SRO  long-  and  short-period 
data  for  7  events  in  this  area  during  1977-79,  whose  locations  and  origin  times  are  given  in 
Table  III -1 .  We  have  concentrated,  in  particular,  upon  the  last  3  events,  given  as  numbers  7, 
8,  and  9  in  the  table. 


TABLE  lll-l 

PRELIMINARY  DETERMINATION  OF  EPICENTER 
LOCATIONS  OF  THE  EVENTS  STUDIED 

— 

— 

Time 

Latitude 

Longitude 

mb 

Event 

Date 

(H:m:s) 

(°N) 

(°N) 

3 

08/29/78 

02:37:06.3 

50.14 

78.16 

5.60 

4 

09/15/78 

02:36:57.3 

49.90 

78.93 

6.00 

5 

1 1/04/78 

05:05:58. 1 

50.02 

79.02 

5.50 

6 

1 1/29/78 

04:33:03.6 

50.00 

78.51 

6.20 

7 

06/23/7? 

02:56:58.6 

49.94 

78.97 

6.20 

8 

07/07/79 

03:46:58.3 

50.06 

79.11 

5.80 

9 

08/04/79 

03:56:58.0 

49.89 

78.96 

6.00 

_ 1 

1.  Short -Period  Data 

The  short -period  P-waveforms  for  all  the  events  studied  were  remarkably  similar  in  both 
amplitude  and  waveform  shape.  Short-period  P-waveforms  observed  at  CHTO  (Chiengmai, 
Thailand),  MAJO  (Matsushiro,  Japan),  ANMO  (Albuquerque,  New  Mexico),  and  KAAO  (Kabul, 
Afghanistan)  are  shown  in  Fig.III-1  for  events  7,  8,  and  9.  Only  very  minor  differences  in 
shape  and  amplitude  can  be  seen,  and  these  similarities  persisted  at  other  stations  and  for 
other  events. 

2.  S-W'aves 

We  also  searched  for  short-period  S-waves,  but  at  all  stations  except  KAAO,  the  detector 

had  switched  off  and  the  data  were  not  recorded.  No  clear  short-period  S-waves  were  observed 
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at  Kabul  tor  any  of  the  events  available.  Long-period  S-waves  were  clearly  observed  at  the 
closer  stations,  and  were  predominantly  SV  motion.  There  was  no  apparent  difference  in  the 
amount  of  SV  generated  from  event-to-event. 

3.  Rayleigh -Wave  Amplitudes 

Large  differences  were  observed  between  the  amplitudes  of  the  Rayleigh  waves  generated 
bv  event  8  and  all  the  others.  Figure  III-2(a-c)  shows  these  amplitude  differences  as  observed 
in  the  time  and  frequency  domains.  The  spectral  amplitudes  have  been  corrected  for  geomet¬ 
rical  spreading  and  attenuation  of  0.75  x  10  4/km,  to  a  common  distance  of  30°;  the  time-domain 
amplitudes  are  directly  measured  from  the  data.  Excluding  event  8,  the  variation  of  time-  or 
frequency-domain  amplitudes  with  azimuth  is  very  similar.  Event  8  has  large  amplitudes  at 
stations  MAJO,  KAAO,  and  BCAO  and  small  amplitudes  at  the  other  stations.  Clearly,  the 
source  mechanism  for  this  event  is  radically  different  from  the  others. 

4.  Rayleigh-Wave  Phases 

In  view  of  the  large  amplitude  differences  observed  for  event  8,  which  appear  to  indicate  a 
lobed  radiation  pattern,  we  have  looked  for  phase  differences  between  event  8  and  the  others. 
Differences  in  event  location  become  significant  here:  for  example,  a  difference  of  15  km  in 
source-receiver  distance  corresponds  to  ~ir/2  phase  shift  at  the  shorter  periods  (t  =  20  sec, 

C  =  3  km/sec,  X  =  60  km).  It  is  thus  important  to  (a)  obtain  the  best  possible  locations  and 
(b)  correct  the  seismograms  to  be  at  a  common  distance  for  each  station.  We  have  used  the 
relative  locations  obtained  for  these  events  as  described  ei-ewhere  in  this  report,  and  cor¬ 
rected  the  spectral  phase  for  each  seismogram  so  that  the  t  it  appears  to  be  at  a  common 
point  chosen  to  be  at  the  approximate  geographical  center  of  the  relocated  events  (50.0°N, 

78.90 °E).  The  furthest  event  from  this  point  (12  km  to  the  southwest)  was  event  9.  We  used 
the  dispersion  for  a  shield  model5  to  adjust  the  phase,  and  the  maximum  phase  delay  added  (or 
subtracted)  was  0.2  circle.  The  equalized  seismograms  were  then  plotted  out  for  all  events 
recorded  at  each  station.  Figure  III-3  shows  seismograms  for  events  7,  8,  and  9  as  recorded 
at  stations  KAAO  and  MAJO,  and  it  can  be  seen  that  those  for  event  8  are  in  both  cases  reversed 
with  respect  to  those  for  events  7  and  9.  This  appears  to  be  a  reversal  only,  occurring  over 
the  entire  frequency  range  of  the  data,  and  no  time  delay  seems  to  be  involved.  Similarly  clear 
reversals  were  observed  at  some  of  the  other  stations,  but  in  some  cases  the  signal  amplitudes 
were  so  reduced  for  event  8  that  phase  changes  were  hard  to  determine.  We  therefore  chose 
to  determine  phase  shifts  with  respect  to  event  9  (the  largest  of  the  events  with  a  "normal" 
radiation  pattern)  in  the  frequency  domain.  At  each  station  the  cross -coirelogram  of  the  seis¬ 
mogram  for  each  event  with  those  for  event  8  was  computed:  its  phase  is  the  phase  difference 
of  the  two  seismograms.  Table  III-2  gives  the  phase  differences  determined.  In  all  cases,  the 
seismograms  for  event  8  v/ere  reversed  (phase  difference  of  ir)  with  respect  to  those  for  event  9. 
In  addition,  the  seismograms  for  the  other  events  were  apparently  similar  in  phase  to  event  9, 
only  in  a  few  cases  did  the  phase  differences  exceed  tt/5. 

5.  Love -Wave  Amplitudes  and  Phases 

Figure  ID -4  shows  the  maximum  Love-wave  amplitudes  observed  on  the  transverse  (after 
rotation)  component  seismograms.  In  all  cases,  the  Love  waves  were  roughly  comparable  in 
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TABLE  II 1-2 

1 

I 

PHASE  DIFFERENCES  (MULTIPLES  OF  it)  MEASURED  FROM  CROSS-CORRELOGRAMS 

OF  EACH  SEISMOGRAM  WITH  THAT  FOR  EVENT  9  (EXCEPT  KONO,  NOT  OPERATING 

FOR  EVENT  9; 

EVENT  4  USED) 

_ 

Station 

Event 

ANMO 

ANTO 

BCAO 

CHTO 

GRFO 

KAAO 

KONO 

MAJO 

SHIO 

3 

+  0.8 

- 

- 

0.0 

- 

+  0.5 

- 

- 

+  0.2 

4 

0.0 

-0.1 

- 

- 

- 

-0.2 

0.0 

-0. 1 

+  0.4 

5 

-0.3 

0.0 

- 

-0.2 

-0.2 

+  0.5 

-0.1 

+  0.6 

0.0 

6 

-0.1 

-0.1 

- 

-0.1 

-0.2 

-0.1 

-0.1 

+  0.2 

0.0 

7 

-0.2 

-0.5 

+  0.2 

-0.2 

- 

+  0.3 

+  0.2 

0.0 

0.0 

8 

+  0.9 

+  0.9 

+  0.8 

- 

+  0.8 

+  0.8 

+  0.7 

+  0.8 

+  0.5 

9 

0* 

0* 

0* 

0* 

0* 

0* 

- 

0* 

0* 

*  By  definition. 

amplitude  with  the  Rayleigh  waves,  and  the  major  difference  again  is  for  event  8  for  which  the 
amplitudes  were  greater  at  all  azimuths  compared  with  the  other  events.  In  addition,  there  is 
an  apparent  radiation  pattern  for  the  Love  waves  which  is  very  similar  in  pattern  to  that  for  the 
Rayleigh  wa/es,  except  that  the  "nodes"  of  the  pattern  are  not  as  well  pronounced. 

We  applied  the  same  equalization  procedure  described  above  for  Rayleigh  waves  to  the  Love- 
wave  seismograms,  and  searched  for  phase  differences  between  the  seismograms  for  different 
events.  Because  of  the  slightly  lower  signal-to-noise  levels  of  the  Love  waves,  the  results 
were  somewhat  less  reliable,  but  no  consistent  phase  differences  were  observed  for  any  of  the 
events.  In  particular,  no  phase  reversal  was  observed  for  event  8  compared  with  the  others. 

6.  Conclusions 

A  careful  analysis  of  SRO  data  recorded  from  7  presumed  explosions  in  a  very  small  area 
(15  X  20  km)  of  the  Eastern  Kazakh  test  site  revealed  large  differences  between  an  event  on 
7  July  1979  and  all  the  others.  These  variations  included  enhanced  Love-wave  generation,  sub¬ 
stantial  azimuthal  variations  in  both  Love-  and  Rayleigh-wave  amplitudes,  and  an  apparent 
phase  reversal  of  the  Rayleigh  waves  at  all  observing  azimuths  and  over  the  entire  period  range 

of  15  to  50  sec  studied.  However,  there  is  no  indication  of  a  time  delay  of  the  surface  waves, 

2  3 

as  observed  by  others.  '  No  appreciable  differences  were  observed  in  the  short-period  P-wave 
data  or  in  the  generation  of  long -period  SV. 

These  differences  clearly  indicate  a  different  mechanism  for  the  event  of  7  July  1979  as 

compared  with  the  5  other  events  nearby.  Possible  explanations  for  these  differences,  which 

2  14 

have  been  previously  observed  for  other  events  both  in  this  and  other  ’  areas,  include  "tec¬ 
tonic  stress  release"  and  spall  closure.  We  note  here  that  neither  seem  particularly  attractive 
explanations:  spall  closure  because  there  are  no  apparent  differences  in  the  short-period 
P-waves,  and  tectonic  stress  release  because  several  of  the  "normal"  evei  cs  are  within  a  few 


kilometers  of  the  anomalous  ones.  Later  in  this  report  we  analyze  the  Rayleigh-wi  ve  radiation 
in  order  to  interpret  these  differences  in  terms  of  the  seismic  moment  tensor. 

R.  G.  North 
M.  A.  Tiberio 

B.  R  BLOC  AT  ION  OF  PR  FSU  MED  EXPLOSIONS  NEAR  THE  NUCLEAR 

TESTING  (.ROUND  IN  EASTERN  KAZAKHSTAN 

There  are  two  principal  motivations  for  this  work,  foremost  of  which  is  a  need  to  normalize 
surface  waves  from  these  events  to  a  common  focal  point  (si  a  Sec.  A  above).  The  other  motiva¬ 
tion  stems  from  a  desire  to  investigate  systematic  errors  in  locations  computed  by  the  usual 
travel-time  analysis.  This  requires  locations,  in  this  case  epicenters,  that  are  refined  in  such 
a  way  as  to  be  substantially  fiee  of  systematic  errors,  the  two  prime  sources  of  which  are  dif¬ 
ferences  between  gross  earth  and  regionalized  travel-time-distance  tables  and  nonrandom  errors 
in  reading  arrival  times. 

Adjustments  toa  common  foc  us  can  be  represented  to  agood  approximation  byexp[ik(x  -“x  )), 
where  k*  x,  and  xc  are  the  wave  number,  location,  and  common-focus  vectors,  respectively. 

To  make  such  adjustments  accurately  requires  locations  that  are  correct  to  one -fifth  to  one- 
tenth  of  the  shortest  wave  length  retained  in  the  spectrum,  which  for  the  digital  SRO  long-period 
data  is  about  15  sec  or  a  wave  length  of  about  50  km  for  fundamental -mode  surface  waves.  Con¬ 
sequently,  locations  must  be  accurate  to  5  to  10  km.  To  measure  accuracy  iather  than  the  pre¬ 
cision  of  the  computation,  there  must  be  at  least  one  reference  event  with  a  known  location. 
Fortunately,  such  an  event  occurred  near  the  eastern  edge  of  this  testing  ground.  It  is  seen  as 
a  crater  on  Landsat  imagery/3  This  crater  is  thought  to  be  the  result  of  an  explosion  on  15  Janu¬ 
ary  1965.  The  crater  lies  at  the  northern  end  of  a  man-made  reservoir  shown  in  Fig.  III-5. 

Accurate  epicenters  were  the  result  of  a  two-step  process.  First,  relative  epicenters 
were  computed  using  arrival  times  of  P  and  PKP  phases  reported  in  the  Preliminary  Determi¬ 
nation  of  Epicenters  (PDEs)  of  the  USGS,  or  times  reported  in  the  ISC  Bulletins.  (We  are  grate¬ 
ful  to  Russ  Needham  for  supplying  us  with  recent  PDEs  prior  to  publication.)  The  relative  loca¬ 
tions  were  computed  by  a  scheme  described  in  the  March  1979  SATS.7  This  scheme  is  well 
suited  to  the  situation  where  distances  between  events  are  much  less  than  the  shortest  ray  path, 
as  is  the  case  for  this  testing  ground.  The  second  step  involves  a  constant  vector  displacement 
wf  about  7  km  toward  the  northwest.  This  displacement  places  the  location  of  the  cratering  event 
(relative  to  that  of  the  master)  at  the  site  of  the  crater.  The  resulting  epicenters  and  standard 
errors  are  listed  in  Table  I1I-3. 

The  master  earthquake  was  chosen  arbitrarily  from  among  the  better  recorded  of  these 
presumed  explosions.  The  event  of  23  November  1976  is  the  master.  The  cratering  event  was 
not  chosen  as  the  master  because  it  was  the  least-well-recorded  of  these  events.  The  number 
of  common  stations  between  the  master  and  the  other  events  ranged  from  32  for  the  cratering 
event  to  136  for  the  11  January  1978  event.  The  closest  stations  were  about  20°  away  near  the 
southern  border  of  the  Soviet  Union.  Stations  internal  to  the  Soviet  Union  do  not  report  events 
within  their  testing  grounds.  The  standard  error  in  the  relative  location  of  the  cratering  event 
has  been  added  to  the  standard  errors  of  the  other  relative  locations  to  yield  the  error  bars 
shown  in  Fig.  Ill- 5.  The  standard  errors  are  in  every  case  less  than  5  km,  which  is  more  than 
adequate  for  the  normalization  of  the  surface  waves  to  a  common  focal  point. 
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TABLE  111-3 
REFINED  EPICENTERS 


Event 

Date 

Latitude 

(°N) 

Longitude 

(°E) 

1 

6/11/78 

49.918  ±0.026 

78.779  ±0.049 

2 

7/5/78 

49.892  ±0.029 

78.  837  ±  0.047 

3 

8/29/78 

50.010  ±0.043 

78.975  ±0.041 

4 

9/15/78 

49.923  ±0.028 

78.842  ±  0.051 

5 

1 1/4/78 

50.027  ±0.039 

78.921  ±0.044 

6 

1 1/29/78 

49.894  ±0.041 

78.898  ±0.057 

7 

6/23/79 

49.894  ±0.032 

78.828  ±0.062 

8 

7/7/79 

50.001  ±0.047 

78. 958  ±0.067 

9 

8/4/79 

49.904  ±0.034 

78.843  ±0.061 

10 

10/29/79 

49.981  ±0.039 

78.945  ±  0.052 

11 

1  1/23/76  (Master) 

50.005  ±0.025 

78.937  ±  0.026 

12 

12/7/76 

49.963  ±0.031 

78.814  ±  0.050 

13 

5/29/77 

49.956  ±0.034 

78.729  ±  0.048 

14 

9/5/77 

50.106  ±0.041 

78.914  ±0.046 

15 

1 1/30/77 

49.975  ±0.032 

78.854  ±0.048 

Cratering 

1/15/65 

49.917 

79.000 

The  arrows  in  the  figure  point  from  the  PDE  location  to  the  refined  location.  As  can  be 
seen,  the  apparent  size  of  the  testing  ground  has  been  reduced  to  a  box  20  km  in  longitude  to 
15  km  in  latitude.  Excluded  from  the  box  is  the  northernmost  event.  It  seems  clear  that,  even 
with  large  sets  of  arrival  times  from  well-recorded  events,  locations  based  on  the  usual  travel¬ 
time  analysis  can  be  significantly  improved  by  some  scheme  that  is  less  sensitive  to  systematic 
errors  in  standard  time-distance  tables  and  in  the  reading  of  arrival  times.  It  should  be  em¬ 
phasized  that  taking  arrival-time  differences  as  was  done  for  the  relative  locations  reported 
here  lessens  the  chance  of  masking  arrival-time  errors  by  a  shift  in  the  computed  location. 

T.  J.  Fitch 
R.  G.  North 

C.  MOMENT-TENSOR  INVERSION  FOR  VERY  SHALLOW  SOURCES 

In  attempting  to  invert  the  surface  waves  generated  by  the  explosion  sources  described  in 
the  previous  section,  serious  problems  emerged  which  are  apparently  inherent  in  moment- 
tensor  analysis  of  sources  at  very  shallow  depths.  In  particular,  the  solutions  for  the  full 
moment  tensor  were  completely  dominated  by  components  Mxz  and  M^.  Numerical  experi¬ 
ments  involving  the  inversion  of  synthetic  seismograms  revealed  the  reasons  for  this. 


-’?=  -.-c— 
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The  far-field  spectra  of  Rayleigh  [FR(h,  w)  ]  and  Love  lF^(h,  u>)]  waves  as  given  by 
Mendiguren  can  be  phrased  as 

FR(h,  u>)  =  SR  {Ar  [  Mxy  sin  2G  -  (M  -  Mxx)  cos  Ee/2] 

+  CRMm  +VM„  +Mx*>/2 
-iBn  [M  sinO  +  M  coso]} 

FL(h,  a>)  =  SL  {Al  [  M  cos  20  +  (M  -  Mxx)  sin  20/2] 

+  i BL  [  -  Mvz  cos  ()  +  Mzx  sin  G  ]} 

v  ere  h  =  source  depth,  w  =  frequency,  ()  =  azimuth,  and  (right-hand  side  notation  that  of 

g 

Mendiguren  ) 


SR  = 

yR(0> 

4CUI4 

f  2 

1  Trkr1 

I*/2  exp[-i(kr 

4>l 

SL  = 

, 

4CUIt  1 

Irtr' 

*/2  exp[-i(kr 

+  7>l 

ar 

R/u. 

y3  (h) 
c 

y4R(h) 
BR  "  pu> 

;  CR 

y2R(h)  +  XkyR(h) 
o;(X  +  2p) 

y/U) 

y^(h) 

II 

J 

< 

c 

• 

bt  =  — — 

L  [ICC 

• 

A  step  function  is  assumed  as  the  source-time  function.  It  can  be  seen  from  the  above  that 
at  any  given  frequency  w,  the  relative  contribution  of  any  given  component  of  the  moment  tensor 

is  governed  by  the  associated  azimuthal  (0)  variation  and  the  values  of  A,  B,  and  C.  The  y.^L 

L  R  R  * 

appearing  in  the  latter  are  the  horizontal  (y.  ,  y .  )  and  vertical  (y?)  displacements,  and  shear 
L  R  R  ii  j 

(y7  ,yA)  and  normal  (y,  )  stresses  obtained  by  solution  of  the  dispersion  equation.  A  necessary 

c  H  ^  L  R  R 

boundary  condition  of  the  solution  is  that  the  stress  components  (y^  .  y2  #  y4 )  vanish  at  the  free 

surface.  Thus,  at  the  surface 

br  “  BL  -  0 

and,  additionally, 

p  -  h.  A  ^  J_  A 
^R  ‘  \  +  2|i  aR  3hR  * 

As  we  go  from  the  surface  waves  downward,  the  stress  components  gradually  increase;  but 
at  very  shallow  depths  h,  such  that  (h/ wavelength)  «  1,  the  above  are  approximately  true. 

o 

We  have  generated  AR,  BR,  and  CR  at  a  variety  of  depths  in  a  shield  model,  and  these  are 
shown  in  Fig.  Ill- 6  over  a  period  range  6  to  50  sec.  At  the  very  shallow  depths  (0.5  and  1.0  km), 
Br  is  much  smaller  than  either  AR  or  CR  except  at  the  shorter  periods  where  surface-wave 
radiation  is  highly  attenuated  for  teleseismic  paths.  Thus,  any  noise  whatsoever  will  be  trans¬ 
lated  into  large  values  of  the  tensor  components  M  ,  M  associated  with  Bn.  In  addition,  since 
p  yz  ** 

the  stress  component  y£^  of  CR  is  very  small,  CR  has  the  same  variation  with  frequency  as  AR, 
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making  it  difficult  to  separate  the  azimuthally  independent  part  of  F^(h,  u>).  This  situation  is 

even  more  accentuated  for  the  Love-wave  variation  of  A.  ,  13.  given  in  Fig.JIl-7.  and  reveals 
10  l,  L, 

the  well-known  lack  of  depth  resolution  provided  by  Love-wave  spectra. 

To  test  the  effects  of  the  above  features  of  A^,  B^,  and  we  have  generated  synthetic 

seismograms  at  a  simulated  network  of  3  stations,  designed  for  good  azimuthal  coverage,  for 


a  source  of  =  1.0  and  all  other  =  0,  corresponding  to  a  vehicle  strike-slip  fault,  at  a 
variety  of  source  depths.  Inversion  was  carried  out  in  the  time  domain  and  the  matrix  decompo¬ 
sition  was  carried  out  in  double  precision.  It  became  apparent  that  even  numerical  noise,  cor¬ 
responding  to  loss  of  precision,  produces  erroneous  results  at  very  shallow  depths.  With  actual 
(rather  than  synthetic)  data,  any  noise  will  automatically  result  in  complete  dominance  of  the 

moment  tensor  by  the  dip-slip  components  of  M  and  M  .  In  order  not  to  exclude  the  pos- 

xz  yz 

sibility  of  such  explosive  sources,  we  did  not  apply  the  constraint  £M..  =  0. 

The  result-  ,  given  in  Tables  III -4  and  III-5,  are  also  simply  summarized  by  Fig.  III-8  which 
shows  the  range  of  eigenvalues  as  a  function  of  depth.  This  can  be  seen  to  decrease  rapidly  from 
2  x  to4  at  0.5  km  to  16  at  10  km  depth.  Constructing  the  inverse  from  the  entire  set  of  eigen¬ 
values  produced  solutions  for  the  seismograms  generated  for  M  ^  0  only  containing  large 

xy 

amounts  of  Mxz  and  M^z  at  0.5  km,  rapidly  decreasing  as  the  source  depth  was  increased. 
Computed  seismograms  for  the  moment  tensor  thus  obtained  were  compared  with  the  originals 
and  could  be  exactly  overlaid.  Seismograms  were  also  computed  for  a  source  corresponding 
to  the  difference  between  the  input  and  output  moment  tensors  (mainly  Mxz,  M^z)  and  these  pro¬ 
duced  reasonable -looking  seismograms  which  were,  however,  a  factor  of  at  least  50  smaller 
than  the  irput  data.  Dropping  the  smaller  eigenvalue  results  in  loss  of  resolution  of  Mzz,  which 
is  traded  with  M  and  M  .  The  exact  solution  (to  4  significant  figures)  is  obtained  at  all  depths 
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only  when  the  3  largest  eigenvalues  are  retained,  but  in  this  case  Mzz  is  very  poorly  resolved, 

and  Mxz,  M^z  essentially  unresolved.  We  conclude  that  for  very  shallow  sources,  only  M^, 

M  ,  and  M  can  be  adequately  resolved.  If  the  constraint  =  0  is  applied,  as  would  be 

\>  yy  zz 

appropriate  for  shallow  earthquakes,  M  equals  -(M  +  M  )  but  little  can  be  done  about  M 

AA  \  y  xz 


and  M  with  the  teles eismic  pass  band. 

yz 


R.  G.  North 
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D.  MOMENT-TENSOR  INVERSION  OF  RAYLEIGH  WAVES 

FROM  EASTERN  KAZAKH  EXPLOSIONS 

In  previous  sections  we  described  (1)  large  differences  in  surface-wave  radiation  between 
presumed  explosions  in  Eastern  Kazakh,  and  (2)  fundamental  problems  in  the  inversion  of  surface - 
wave  data  from  very  shallow  sources  to  determine  the  full  moment  tensor.  We  attempt  here, 
subject  to  the  limitations  posed  by  the  latter  problem,  to  invert  the  Rayleigh  waves  in  order  to 
characterize  the  nature  of  the  source  in  terms  of  the  resolvable  components  of  the  moment  tensor. 

For  an  isotropic  explosive  source,  the  moment  tensor  should  be  characterized  by  equal 
diagonal  elements  (M^  =  M^  =  Mzz),  and  the  off-diagonal  elements  should  be  zero.  The  fact 
that  Love  waves  are  virtually  always  generated  by  explosions  indicates  that  such  a  simple  de¬ 
scription  of  the  source  is  inadequate.  An  isotropic  explosive  source  should  also  produce  iden¬ 
tical  seismograms  at  all  azimuths;  this  is  difficult  to  determine  from  data  due  to  our  lack  of 
knowledge  of  propagation  corrections,  but  some  explosions  do  produce  what  seem  to  be  azi- 
muthally  dependent  radiation  patterns. 
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TABLE  111-4 

SOLUTIONS  TO  SYNTHETIC  M  SEISMOGRAMS* 

xy 


Source  Depth 
(km) 

M 

XX 

M 

xy 

M 

xz 

M 

yy 

M 

y* 

M 

zz 

0.5 

0.15 

1.0 

0.56 

0.15 

0.35 

0.49 

1.0 

0.027 

1 .0 

0.27 

0.027 

0.17 

0.094 

1.5 

0.009 

1.0 

0.17 

0.008 

0.11 

0.032 

2.0 

0.003 

1.0 

0.12 

0.075 

0.075 

0.013 

3.0 

0.0 

1.0 

0.074 

0.046 

0.046 

0.003 

5.0 

0.0 

1.0 

0.036 

0.022 

0.022 

0.003 

10.0 

0.0 

1.0 

0.007 

0.006 

0.006 

0.001 

*  Solution  vector  obtained  by  inversion  of  MXy  =1.0,  other 

Mjj  =0,  as  a  function  of  depth,  when  all  six  eigenvalues 
are  used  to  construct  the  inverse. 

TABLE  II 1-5 

RESOLUTION  MATRIX* 

M  M 

xx  xy 

M  M 

XZ  yy 

M 

y* 

M 

zz 

1.0  0.0 

0.0  -0.078 

0.0 

-0.26 

M 

XX 

1.0 

0.0  0.0 

0.0 

0.0 

M 

xy 

0.0  0.0 

0.0 

0.0 

M 

xz 

0.921 

0.0 

-0.26 

M 

yy 

0.0 

0.0 

M 

y* 

0.16 

M 

yy 

*  Resolution  matrix  obtained  for  0.5-km  depth 
3  eigenvalues  only  are  retained.  Note  total  l< 
MXZ/MyZ  and  very  poor  resolution  for  Mzz,  w 
trade-off  with  M  and  Mw. 
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source  when  largest 
ack  of  resolution  for 
hich  has  considerable 

A  major  problem  in  moment-tensor  inversion,  or  anv  other  attempt  to  characterize  the 
source,  is  separation  of  path  and  source  effects.  These  problems  are  particularly  severe  when 
short-period  (up  to  50  see)  surface  waves  are  used  since  these  sample  mainly  the  uppermost 
mantle  and  crust,  where  major  lateral  inhomogeneities  exist.  At  shorter  periods  (~20  sec), 
phase  and  group  velocity  variations  of  up  to  25  percent  are  quite  common,  and  attempts  to  de¬ 
termine  attenuation  for  all  but  the  simplest  paths  frequently  produce  such  disastrous  features 

1 1 

as  negative  QJ  It  has  been  pointed  out  that  correct  knowledge  of  the  phase  velocity,  to  better 

than  0.5  percent,  is  necessary  to  determine  source  phase  to  sufficient  accuracy  for  inversion. 

To  minimize  the  effects  of  poor  knowledge  of  path  characteristics,  we  have  used  a  reference 

1 1 

event  approach  somewhat  similar  to  that  of  Patton. 

Ideally,  the  characteristics  of  propagation  paths  from  a  particular  source  region  should  be 

determined  from  the  surface  waves  generated  by  a  large  event  of  known  moment  tensor  in  that 

region.  The  path  corrections  thus  obtained  can  then  be  applied  to  the  radiation  from  smaller 

events  to  determine  their  source  moment  tensor.  The  moment  tensor  of  this  calibration  event 

can,  for  example,  be  estimated  from  the  fault-plane  solution  for  the  source  if  it  is  large  enough 

to  generate  the  long-period  P-waves  required.  This  is  not  feasible  in  the  present  case,  as  the 

moderate-size  explosions  studied  were  very  poor  generators  of  long-period  P. 

We  have  thus  chosen  to  invert  for  the  four  resolvable  components  (M  ,  M  ,  M  ,  M  )  of 

xx  yy  zz  xy 

the  moment  tensor  in  the  following  two  ways: 

(1)  Amplitudes  only,  and 

(2)  Differences  between  the  seismograms  for  each  event  and  those  for 
event  9,  which  we  have  selected  as  a  calibration  event.  While  we 
cannot  determine  the  absolute  values  of  phase  delay,  we  earlier 
determined  the  phase  differences  at  each  station  between  each 
event  and  event  9. 

1.  Data -Reduction  Procedure 

The  first  step  here  involved  correcting  some  obvious  errors  in  the  data.  The  timing  at 
KAAO  (Kabul)  was  clearly  off  by  exactly  1  min.  for  late  1978  and  early  1979:  this  was  corrected. 
No  anti-alias  filter  had  been  applied  to  the  long-period  vertical  recordings  at  GRFO:  the  data 
were  low-pass  filtered  to  make  all  seismograms  comparable  in  instrument  response.  Some  of 
the  stations  are  SROs  (CHTO,  SHIO,  ANMO,  ANTO,  TATO,  GRFO,  BCAO):  the  others  (MAJO, 
KAAO,  KONO)  are  ASROs.  For  uniformity,  the  seismograms  at  the  ASROs  were  corrected  to 
appear  as  if  they  were  recorded  on  SRO  instruments:  this  involved  corrections  in  the  frequency 
domain  for  the  small  differences  in  amplitude  and  phase  response  between  the  two  types  of 
instruments. 

For  some  of  the  events,  there  was  some  interference  from  other  sources:  to  reduce  this 
and  improve  signal-to-noise  ratios,  we  determined  group-velocity  curves  for  each  path  from 
the  event  with  the  largest  and  simplest  signals  and  used  these  to  design  a  time  variable  filter. 
Particular  care  was  taken  to  ensure  that  phase  differences  between  events  were  preserved. 

The  seismograms  had  all  previously  been  equalized,  as  described  in  Sec.  A  above,  to  a  common 
focal  point  of  50.0° N,  78.9°E. 

We  have  chosen  to  carry  out  the  inversion  in  the  time  domain.  To  construct  seismograms 
containing  amplitude  information  only,  we  have  set  the  phase  delay  to  be  simply  (^  +  j  - 1-)  where 
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—  is  the  propagation  delay  (A  =  distance  (km);  c  =  phase  velocity)  and  the  and  are  re¬ 
quired  to  account  for,  respectively,  the  assumed  step  source-time  function  and  the  effects  of 
superposition  in  the  frequency  domain  to  construct  the  seismogram.  For  the  difference  in  seis¬ 
mograms  with  respect  to  those  for  event  9,  we  have  (a)  added  an  additional  phase  delay  which 
is  the  phase  difference  between  the  seismogram  and  that  for  the  same  path  from  event  9,  and 
(b)  constructed  seismograms  which  are  the  difference  between  these  and  those  for  event  9. 

2.  Inversion 

As  noted  earlier,  we  carry  out  inversion  in  the  time  domain.  This  involves,  for  each  set 
of  seismograms  for  a  particular  event,  calculating  6  synthetic  seismograms  -one  for  each 
component  of  the  moment  tensor  -  and  then  inverting  to  determine  the  optimum  linear  combina¬ 
tion  of  these  which  fit  the  data  in  a  least -squares  sense.  The  linear  relation  of  the  moment 
tensor  to  the  seismogram  means  that  inversion  of  difference  seismograms  will  yield  the  differ¬ 
ence  moment  tensor.  The  synthetic  seismograms  are  computed  using  the  solution  at  1  km  depth 
in  a  shield  model.  The  dispersion  for  the  same  model  is  applied  as  a  path  correction  to  both 
the  data  and  the  synthetics:  the  choice  of  model  here  is  not  important,  since  it  is  canceled  in 
the  solution  procedure.  Since  we  know  that  the  sources  are  by  definition  very  shallow,  we  do 
not  need  to  include  source  depth  in  the  inversion  procedure.  The  depth  of  1  km  used  may  not 
be  correct,  but  from  Sec.  C  above  we  are  aware  that  the  excitation  parameters  and  do 
not  change  appreciably  over  the  range  0  to  1  km  depth  and  that  the  most  important  feature  of  the 
solution  is  that  is  very  small,  giving  extremely  poor  resolution  for  tensor  elements  Mxz 

and  M  .  We  are  thus  able  to  solve  only  for  the  remaining  components  M  ,  M  ,  M  ,  and 
yz  xx  yy  zz 

Mxv>  with  the  additional  complications  of  trade-off  between  the  diagonal  elements  Mxx«  M  , 
and  M  .  The  interpretation  of  such  a  solution  is,  of  course,  extremely  limited. 

3.  Results  with  Amplitudes  Only 

The  6  synthetic  seismograms  corresponding  to  individual  moment-tensor  components  are 
computed  with  the  source  phase  set  to  zero.  The  attenuation  coefficients  used  to  correct  ampli¬ 
tudes  are  given  in  Table  III-6.  Constructing  the  solution  from  the  3  largest  eigenvalues  (those 
rejected  correspond  to  unresolved  Mxz>  MyZ,  and  poor  resolution  of  Mzz),  the  solutions  are 
shown  in  Table  III-7  for  events  4  through  9.  The  difference  between  event  8  and  the  others  is 
immediately  apparent  from  the  solution  -  M  is  much  larger,  and  reversed  in  sign,  for  this 

xy 

event  compared  with  the  others.  If  we  were  able  to  resolve,  and  diagonalize  the  full  moment 
tensor,  then  the  explosive  or  monopole  component  would  be  given  by  the  diagonal  elements,  and 
the  multipole  component  by  the  other  elements.  The  ratio  of  M  to  the  sum  of  the  diagonal  ele- 
ments  is  also  given  in  Table  III-7:  it  is  much  larger  and  reversed  in  sign  for  event  8  compared 
with  the  others.  In  the  following  contributions  to  this  SATS,  the  problem  of  separating  the  multi- 
pole  components  of  these  incomplete  moment  tensors  is  discussed  and  possible  solutions  to  this 
problem  are  proposed.  The  difference  seismograms,  with  phase  differences  included,  are  of 
particular  use  in  this. 

4.  Comparison  of  Predicted  and  Observed  Radiation  Patterns 

Using  the  4  resolved  components  of  the  moment  tensor,  we  may  compute  radiation  patterns 
for  comparison  with  the  variation  of  amplitudes  with  azimuth  described  in  Sec.  A  above.  In  the 
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absence  of  components  and  Myz  we  may  ^determine  only  the  real  part  of  the  far-field  spec 
trum.  given  by 

[  <Myy  -  , 

FR(u)  -  SR(u>)  Ar(w)  Im  sin  20  — — —  cos  20 

-  I  M«  *  (Myv  1  Mxx»/2] 
iur  Rayleigh  waves,  and 


<Mw~M  xx) 

Fl(w)  =  SL(co)  Al(w)  Mxy  cos  20  +  — ZZ-j - —  sin  20 


where  the  nomenclature  is  that  of  Sec.  C  above  and  we  have  assumed  CR  =  -l/3  AR  for  a  very 
shallow  source.  Except  for  the  scaling  introduced  by  SR(u>),  AR(u>)  and  S^(u;),  A^(u’)  the  radia¬ 
tion  patterns  obtained  are  independent  of  frequency. 

The  spectral  amplitudes  calculated  at  azimuths  corresponding  to  the  stations  used  are  shown 
in  Fig.  III-9(a-b)  for  events  7,  8,  and  9,  and  should  be  compared  with  the  observed  variation  into 
azimuth  shown  in  Figs,  in-2  and  III-4.  We  have  inverted  only  the  Rayleigh  wavo^’  to  determine 
the  moment-tensor  components.  The  match  for  the  Rayleigh  waves  is  good:  the  similarity  of 
events  7  and  9  is  clear,  as  is  the  substantial  difference  between  these  two  and  event  8. 

The  major  Love-wave  differences  observed  between  events  7,  e.  and  9  are  (a)  the  larger 
amplitude  at  all  azimuths  for  event  8,  and  (b)  smaller  differences  between  events  7  and  9.  These 
features  are  well  predicted  by  the  moment  tensors,  as  shown  in  Fig.  in-9(b).  The  large  ampli¬ 
tudes  observed  at  KAAO  for  events  8  and  7  are  not,  however,  predicted.  Nevertheless,  it  is 
encouraging  that  the  moment  tensor  determined  only  from  Rayleigh  waves  predicts  the  general 
features  of  the  Love -wave  data.  R  ^  North 

t!  J.  Fitch 


E.  TOWARD  AN  UNDERSTANDING  OF  THE  AZIMUTHALLY  DEPENDENT 

SURFACE-WAVE  RADIATION  FROM  UNDERGROUND  EXPLOSIONS 

From  a  previous  discussion  by  North  and  Fitch  (Sec.  C  above),  a  useful  point-source  repre¬ 
sentation  for  underground  explosions,  if  it  is  to  be  couched  in  terms  of  the  moment  tensor,  must 
be  constrained  if  the  data  are  fundamental  mode  surface  waves.  Solutions  for  the  full  moment 
tensor  are  not  realizable  because  of  the  poor  resolution  of  the  Mzz  and  the  dip-slip  components 
Mxz  and  Myz  from  sources  within  a  few  kilometers  of  the  surface.  At  this  writing,  a  variety 
of  constraints  are  being  tested,  with  only  tentative  indications  of  what  may  prove  to  be  useful 
constraints. 

Trial-and-error  fits  to  a  ’’faultless"  moment  tensor  using  only  amplitude  data  from  Eastern 
Kazakhstan  (see  Sec.  A  above)  suggest  that  there  is  not  sufficient  information  in  these  data  to 
resolve  this  moment  tensor  even  if  it  is  constrained  to  maximize  the  monopole  component.  The 
"faultless"  moment  tensor  in  diagonalized  form  can  be  written 
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a  ami  l)  pertain  to  the  monopole  and  compensated  linear  vector  dipole  (CLVI.))  components,  re¬ 
spectively.  11’  ()  is  the  polar  angle  of  the  major  axis  of  this  tensor  and  is  the  azimuth  of  this 
axis,  the  tensor  components  can  be  written 

■>  ?  2  2  2 
Mxj,  A(cos“0cos  O  4  sin  0)  4  0  cos  0  sin  O 

1/2  sin20(X  cos2o  -  X  +  0  sin2Oj 

2  2  2  2  2 
M  \(sin  </>  cos  O  +  cos  0)  +0  sin  '  0  sin  0 

Mxz  =  1/2  sin  20  cos  0(0  -X) 

M  =  1/2  sin  2o  sin  0(0  -  X) 
yz 

M  =  X  sin20  +  0  cos20 
zz 


Because  this  source  representation  is  axially  symmetric  when  the  M  and  M  terms  are  not 

yz 

resolvable,  there  is  a  180°  ambiguity  in  the  azimuth  of  the  major  axis. 

If  we  take  M  ,  M  .  and  M  as  the  resolved  components  of  the  full  moment  tensor  (see 
yx  xy  yy 

discussion  in  Sec.  D  above),  then  the 
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which  leaves  us  with  the  problem  of  solving  for  X  and  0  given  any  two  of  the  knowns  and  a  trial 
value  of  0.  For  each  trial  value  of  0  in  the  range  0°  to  180°,  the  rms  of  the  fit  is  computed. 
The  results  of  the  Eastern  Kazakhstan  events  revealed  a  very  broad  minimum  in  the  rms  values 
centered  at  a  0  of  90°  which  corresponds  to  a  horizontal  major  axis  and  a  maximum  in  the  mono¬ 
pole  component.  It  seems  unlikely  that  a  horizontal  major  axis  for  the  CLVD  can  represent  the 
actual  physiol  mechanism  for  each  of  these  underground  explosions. 

The  inclusion  of  phase  information  in  the  data  would  substantially  increase  the  quality  of  the 
data;  however,  this  can  not  be  done  at  periods  of  less  than  50  sec  unless  the  propagation  effect 
is  removed.  To  remove  this  effect  a  reference  event  (No.  9  in  Table  m-7)  was  chosen  and  spec¬ 
tral  phase  and  amplitude  differences  were  used  in  inversions  for  the  difference  moment  tensor. 
Since  events  7,  8,  and  9  from  Eastern  Kazakhstan  have  nearly  the  same  body-wave  magnitude, 
the  difference  seismograms  should  be  largely  independent  of  the  monopole  or  simple  explosive 
component  of  the  source  functions.  However,  the  interpretation  of  the  difference  moment  ten¬ 
sors  is  not  simple  because  the  reference  event,  whichever  one  might  be  chosen,  does  not  radiate 
as  a  simple  monopole.  For  example,  there  is  always  an  appreciable  Love  wave  recorded. 

The  difference  moment  tensor  at  a  given  station  is  proportional  to 


-A 


R 


-iQ 


(w)e 


R 

s 


where  A  and  Q  refer  to  amplitude  and  phase,  and  the  subscript  s  and  superscript  R  refer  to 

source  and  reference  event,  respectively.  6Qs  is  the  phase  difference  between  another  event 

and  the  reference  event  at  a  particular  station.  Unfortunately,  the  phase  of  the  reference  event 
R 

by  itself  Q"  has  been  estimated  from  the  M  ,  M  ,  and  M  components  of  the  moment-tensor 
s  xx  yy  ,.y 
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solution  to  the  amplitude  data  from  the  reference  event.  It  is  assumed  here  that  seismograms 
from  all  events  have  been  reduced  to  a  common  focus. 

Table  III- 8  summarizes  the  solutions  lor  the  moment-tensor  differences  of  events  7  and  8 
with  respect  to  event  9.  Notice  that  the  zero  trace  constraint,  which  is  linear,  has  been  used. 
As  in  the  previously  discussed  inversions  using  only  amplitudes,  the  eigenvalues  for  the  dip- 
slip  components  AM  and  AM  are  more  than  one  order  of  magnitude  less  than  those  for  the 
other  three  components.  The  dip-slip  components  are  retained  because  they  don't  trade  off 
with  the  others;  consequently,  you  can  choose  to  believe  them  or  not  without  affecting  the  in¬ 
terpretation  of  the  three  well-resolved  components. 

A  comparison  of  the  AM  ,  AM  ,  and  AM  with  the  corresponding  components  in  Table  III— 7 

xx  jj  xj 

shows  that  inversions  with  and  without  phase  are  fundamentally  different.  With  that  revelation, 
it  becomes  imperative  to  compute  the  phase  of  the  reference  event  in  a  less  arbitrary  manner. 

A  scheme  involving  a  joint  inversion  of  all  the  available  data  from  this  testing  ground  is  now 

being  pursued.  T.  J.  Fitch  A.  M.  Dziewonski 

R.  G.  North  J.  H.  Woodhouse 

F.  THE  SECOND  MOMENT  TENSOR  AND  THE  LOCATION 
OF  SEISMIC  SOURCES 

Specific  formulas  are  given  for  the  calculation  of  theoretical  seismograms  due  to  a  source 
oi  finite  spatial  and  temporal  extent,  including  the  terms  involving  the  second  spatial  and  tem¬ 
poral  moments.  Using  these  results,  expressions  are  derived  for  the  partial  derivatives  of 
synthetic  seismograms  with  respect  to  source  location  and  time  of  occurrence,  which  may  be 
used  for  the  location  of  earthquakes  using  complete  waveform  data.  The  difference  between 
locations  and  times  so  obtained  from  long-period  data,  and  those  obtained  conventionally  from 
short-period  data  can,  if  sufficiently  well  determined,  give  information  about  the  spatial  and 
temporal  extent  of  the  source  process. 

Following  Gilbert,^  Gilbert  and  Dziewonski,^  and  Backus  and  Mulcahy,^  lie  seismic  dis¬ 
placement  of  a  spherically  symmetric  Earth  model  due  to  any  seismic  source  may  be  written 


s(x,  t)  =  Yj  ak(t)  £k(x) 


(HI-1) 
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with 
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S,  L 


1  -  e 


COS  a'k(t  ~  t') 


e^( x f ) : F (x1,  t')  dt'  d  x' 


(III-2) 


The  notation  used  here  and  later  will  follow,  for  the  most  part,  that  of  Gilbert  and  Dziewonski.1 3 
The  vector  fields  s  k  are  a  complete  set  of  eigenfunctions  with  eigenfrequencies  Wj.,  and  the 
label  k  represents  the  four  parameters  (q,  I,  m,  n):  mode  type  (toroidal  or  spheroidal),  an¬ 
gular  order,  azimuthal  order,  and  radial  order,  respectively.  The  tensor  fields  e,  are  the 
elastic  strain  in  the  km  mode,  and  the  overbar  denotes  complex  conjugation.  The  time  deriva¬ 
tive  of  the  stress  glut  tensor  £(x,  t)  completely  characterizes  the  source,  and  the  spatial  inte¬ 
gration  in  Eq. (IH  — 2 )  is  over  the  source  region  V  ,  namely  that  region  in  which  F(x,  t)  is  different 
from  zero.  The  small  parameter  «k  characterizes  the  attenuation  of  the  kin  mode  and  is  obtain¬ 
able  from  perturbation  theory  for  specific  models  of  attenuation  in  the  Earth;  terms  of  fir -t  and 
higher  orders  in  a,  have  been  neglected  in  Eq.  (III-2)  and  will  be  neglected  elsewhere. 

Backus  and  Mulcahy1*  have  shown  that  for  modes  with  wavelengths  long  in  comparison  with 
the  source  dimension,  and  with  periods  long  ir  comparison  with  source  duration,  Eq.  (III-2)  may 

be  usefully  expanded  in  terms  of  the  low-order  spatial  and  temporal  moments  of  the  glut  rate  dis- 

f  -ak(t-f) 

tribution.  Sp  ecifically,  retaining  the  first  two  terms  in  a  Taylor  expansion  of  g,k(x')  [1  -  e  x 

cos  a>k(t  -  t')|  about  a  fiducial  location  =  x,c  in  the  source  region,  and  a  fiducial  time  t1  =  tQ 
close  to  the  origin  time,  Eq.  (Ill -2)  becomes 


ak(t)  =  1  -  e 


.1 


s(k) 


cos  cot(t  -  t0)j  [M..e.\  '(xn)  +A..ne..  Jxn)J 


*ij  ~ij  "-O'  ■  ijp  ij,  pv~0' 

(IB-3) 

at  times  after  the  source  has  ceased  to  act,  where  the  following  glut  moments  have  been  defined: 


sinc^t  ~tQ) 


M..  =  f  \  f..(x,  t)d3xdt 

*'y  J  -oo 


(III -4) 

(in-5) 

(III -6) 


M  is  "the”  moment  tensor  (the  sign  convention  used  here  is  opposite  to  that  of  Gilbert  and  Dzie- 
~  13 

wonski  )  and  A  and  H  describe,  in  the  lowest  approximation,  the  spatial  and  temporal  distribu¬ 
tions  of  the  source  around  the  fiducial  location  xQ  and  the  fiducial  time  tQ.  Following  Backus,15 
we  define  the  centroid  location  and  the  centroid  time  t£  as  the  values  of  Xq,  tQ  for  which  A 
and  H  are  minimum,  in  the  sense  that 


AijpAijp'  HijHij 
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Utain  their  smallest  values.  It  is  readily  shown  that 


A..  M.. 

x  -  v  f  UP  ,-U 

\v  °P  Mp8Mr8 


(111-7) 


V  =  ‘o  + 


The  momcn  tensors  A  and  11  relative  to  this  centroid  will  be  denoted  by 


Aijp  *  AiJP,5c> 
nij  .  H.j(tc) 


and,  thus, 


Aijp(~o)  =  Aijp-(X0p  -V 

VV-flll-(t0-‘c,MiJ  • 


The  tensors  A,  fl,  which  satisfy 


M.  .A . .  =  0 

i]  1JP 

M.  .H.  .  =  0 
i]  i] 


have  the  advantage  of  depending  only  upon  the  source,  and  not  on  the  choice  of  xQ,  tQ. 
Equation  (III-3)  may  ni  be  written 


ak(l)=  1-e 


cos  wk(t  -  tQ) 


x  {My  [S^’tJSo)  +  <*cp  -  V  4>~°)]  +  XMP  5ijkp(~0,} 


-"ke 


sina>k(t  -  t0)  [Mj.ey5*  (x0)(tc  lo*  +  f^ijeij 


(IIi-8) 


(III-9) 


(111-10) 


(II1-11) 


(III— 1 2) 


(111-13) 


( III -14) 


(III— 1  3) 


It  is  clearly  feasible  and  desirable  to  obtain  constraints  upon  A.H  from  long-period  data  for  par¬ 
ticular  earthquakes;  if  it  is  assumed,  however,  that  these  contributions  are  small,  Eq.  UU-1 
leads  to  a  simpler  inverse  problem  for  the  centroid  location  x,  and  the  centroid  time  tc-  e 


ak(t)  =  aok(t)  +  (xcp  -xQp)  5J-  aok(t)  +  dc  t0)  ^  a„k(t) 


(III-16) 


where 


a0k(x't}  =  1  "e 


coswk(t-tQ)  ek(x  q):  M 


(HI-17) 


In  order  to  calculate  Eqs.  (III-3),  (111-15),  or  (III-16),  it  is  necessary  to  obtain  expressions 
for  a.  (t)  in  terms  of  the  precomputed  scalar  eigenfunctions  for  a  particular  Earth  r.iooel.  We 

K  j  o 


then  have 


im</> 


r(x,t)=  Xak!t)Uk(r)X(m(0)e 
k 
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(III-18) 


s()(x.  t)  =  £  ak(t)  l  vk(r)  Xfm<0)  4  im  ct-oco  Wk(r;  xfm(0)|  e1™0 
k 

s0(x>t)  r  Zak(t)  tim  cosec0  Vk(r)  X|xn(0)  -Wk(r)  Xf'm(0)]  eim0 
k 

where  the  notation  follows  that  of  Gilbert  and  Dziewonski.1  3 

We  denote  by  f.  the  six  independent  elements  of  the  moment  tensor:13 


=  Mrr  • 

f2  *  Mee  • 

*-*s 

II 

s 

■©- 

■e- 

f4 

=  MrO  • 

f5  =  Mr*  • 

f6  «  Me* 

and  by  the  corresponding  elements  of  strain  at  the  source  (the  mode  index  k  is  suppressed) 
fl  =  err  *  €2  =  e00  *  €3  =  e00 

€4  =  2ere  •  €5  =  2er  0  '  e6  =  2e00  *  (HI-19) 

Similarly,  we  shall  denote  by  X^,  ei  the  18  independent  elements  of  A,Ve: 

X1  =  Arrr  '  X2  ~  A00r  ’  X3  =  A00r  '  X4  =Ar0r  ’  X5  =  A.0r  ' 


X6  =  A00r  ’ 

X7  =  Arr0  ' 

X8  =  A000  '  X9  =  A 0 00  ’ 

X10  =  Ar00 

’  X11  =  Ar00 

*  X12  =  A000  '  X13  =  Arr0  ' 

X14  =  A000 

'  X1  5  =  A  0  0  0 

’  X16  =  Ar00  '  X17  =  Ar 00 

X18  =  A000  *  (HI-20) 
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TABLE  1 1 1-9 

PARAMETERS  FOR  THE  EXCITATION  DUE  TO  FIRST-  AND  SECOND-MOMENT  TENSORS 


Specific  forms  for  t~.  are  given  by  Gilbert  anti  Dziewonski  J  and  are  repeated  in  Table  III-9. 
Forms  for  e^  also  are  listed  in  Table  1II-9,  and  the  notation  is  given  in  Table  III- 1 0.  The  par 
tial  derivatives  with  respect  to  source  location  and  time  of  occurrence  in  Fq.  (III-i6)  are 


-<*k(t-t0)  l/v  - 

1  -e  cosa>k(t  -t0)  (  7  Te. 


1  0aOk  [.  ‘Q!k(t"t0)  .  ]  (  y  -  \ 

tt  ier=  1  ~e  cosck(t-t0)  2,  tie.+6 


TABLE  111-10 

NOTATIONS  USED  IN  TABLE  111-9 


Spheroidal  Modes 


Toroidal  Modes 


F  =  r  (2U  -1(1  +  1 )  V] 
E$  =  V  -  r"1  (V  -  U) 


ET  =  W  -  r  'W 


Pl  =  r’1(U  -  V  -  r"1  U  +  r"1  V) 
p2  =  r‘2(U  -  2V) 

p3  =  r"1  [U  -  r-1U  +  1/2 r-1 1(1  +  1)(2V  -  U)] 
P4  =  r  V  -  r"1  V) 

p5  =  r"1  (1/21(1  +  1)(V  -  r-1V)  -  U  +  r^U] 

p6  =  r"1  [}/2rm]  i{/t  +  1)V  -  V  -  r"1U] 

p?  =  -r”2  { 1/4  [1(1  +  1)  -2]  V  +  rF} 


=  — r’V  -  r']W) 
q2  =  -2r"2W 


q5  =  1/2  r”1 1  (I  +  1)(W  -  r“V 
q6  =  r”1  (1/2  r-1I(jf  +  1)W  -  W] 
q?  =  -l/4r"2  [1(1  +  1)-2]W 


U,  V,  W  to  be  evaluated  at  radius  r  =  | 


k  -  1/21+1  .  (I  +  n)! 
n  -n  \  4ir  (I  -  n)! 


+  n)!V/2 


k0  “  (^)'/2  !  kn  =  1/2  Ki  -  n)(i  -  n  ♦  1)],/2  kn., 


where  (roi()Q'  ^ 0 ^  arc  ^ie  sPh°rical  coordinates  of  the  source  location  Xq. 


Discussion 

Detailed  formulas  have  been  given  for  the  formulation  of  the  linear  inverse  problem  for 
moments  up  to  the  second,  and/or  for  the  linearized  inverse  problem  of  source  location  using 
complete  waveform  data.  The  problem  of  simultaneously  determining  the  source  moment  ten¬ 
sor  M  and  the  centroid  location  and  time  is  not  linearizable  and  must  proceed  by  iterations  (see 
Sec.G  below)  in  which  the  moment  tensor  and  source  location  are  alternately  fixed,  and  where 
the  process  is  found  to  be  rapidly  convergent. 

It  may  be  remarked  that  it  is  only  in  the  limit  of  long  periods  and  wavelengths  that  the  true 
centroid  location  and  time  will,  even  theoretically,  be  obtained.  In  general,  for  a  finite  number 
of  records  and  a  finite  frequency  band,  a  location  will  result  which  depends  upon  frequency  and 
station  distribution.  Nevertheless,  it  is  to  be  expected  that  for  relatively  small  sources  the 
answer  will  be  close  to  the  centroid,  and  in  any  case  will  give  a  valuable  measure  of  the  size 
and  duration  of  the  source  when  compared  with  conventionally  determined  locations  and  origin 
times,  provided,  of  course,  that  these  are  determined  with  sufficient  accuracy.  If  the  process 
could  be  carried  out  in  a  number  of  different  frequency  bands,  the  resulting  locations  and  times, 
as  a  function  of  frequency,  would  be  very  interesting  data  for  the  study  of  sources. 

J.  H.  Woodhouse 


G.  INVERSION  OF  SRO  AND  ASRO  SEISMOGRAMS  FOR  THE  MECHANISM 

AND  LOCATION  OF  THE  SOURCE 

This  report  describes  application  of  a  novel  procedure:  simultaneous  inversion  of  the  wave¬ 
form  data  for  the  equivalent  system  of  forces  at  the  source  and  the  hypocentral  parameters.  Usu 
ally,  location  of  an  earthquake  and  estimation  of  source  mechanism  are  performed  separately. 
Inversion  for  the  moment  tensor  is  made  assuming  that  the  hypocentral  parameters  are  known. 
This  may  not  always  be  a  correct  assumption;  for  sources  of  finite  size,  the  best  point-source 
location  need  not  be  coincident  with  the  hypocenter  determined  from  the  arrival  times  of  short- 
period  body  waves. 

The  approach  presented  here  stems  from  earlier  work  of  Dziewonski,1^  who  outlined  an  ef¬ 
ficient  procedure  for  evaluation  of  the  elements  of  the  moment  tensor  by  inversion  of  entire 
wavetrains  consisting  of  long-period  body  waves.  A  very  significant  development  has  been  deri¬ 
vation  by  Woodhouse  (see  Sec.  F  above)  of  the  coefficients  for  perturbation  of  excitation  param¬ 
eters  due  to  changes  in  the  hypocentral  parameters  of  the  source  and  in  the  origin  time. 

We  feel  that  our  approach  has  significant  potential  for  research  in  seismology.  In  addition 
to  providing  very  stable,  routine  determinations  of  moment  tensor  for  earthquakes  in  a  range 
of  body-wave  magnitudes  from  5.5  to  7.0,  differences  between  the  standard  hypocentral  locations 
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and  locations  of  the  centroids  of  the  stress  glut  could  lead  to  reliable  estimates  of  temporal  and 
spatial  dimension  of  the  source  region. 

It  is  also  possible  that  the  standard  locations  are  in  error;  we  are  poorly  protected  against 
absolute  errors  in  hypocentral  location  if  lateral  heterogeneity  is  not  known.  It  may  be  expected 
that,  by  using  signals  of  wavelengths  more  than  an  order-of-magnitude  greater,  this  effect  could 
be  detected. 

Another  aspect  of  the  proposed  method  is  its  high  resolution  for  the  source  depth.  Because 
we  are  analyzing  entire  wavetrains  of  body  waves,  the  energy  radiated  upward  and  reflected 
from  the  surface  is,  on  average,  equal  to  the  downward  radiated  energy:  thus,  we  include  the 
complete  suite  of  the  "depth  phases"  in  our  analysis. 

The  process  described  here  requires  determination  of  a  starting  solution  for  the  moment 
tensor  using  hypocentral  parameters  furnished  by  the  National  Earthquake  Information  Service, 
for  example.  This  part  of  the  procedure  is  directly  related  to  the  approach  of  Gilbert  and 
Dziewonski.13 

6 

uR(A,  0,  t)  =  £  A.°  (A.  0,  rg,  t)  •  M.°  (HI-23) 

i=l 

where  u,  is  the  kth  record  in  the  set  of  seismograms  to  be  inverted.  A..  (A,  0,  r  ,  t)  is  the  syn- 
thetic  seismogram  associated  with  the  l  n  component  of  the  unit  moment  tensor,  and  IVh  are  the 
elements  of  the  moment  tensor  (indices  from  1  through  6  are  assigned  to  Mrr,  MQ0,  M^,  M^, 
Mr<£,  and  M respectively).  A.^  are  obtained  by  summation  of  normal  modes;  the  superscript 
denotes  the  zeroth  (starting)  iteration. 

In  practice,  we  use  the  data  from  SRO  and  ASRO  stations.  Theoretical  seismograms  are 

obtained  using  the  normal-mode  catalog  computed  by  Buland17  for  an  Earth  model  1066B  of 
13 

Gilbert  and  Dziewonski.  Since  this  catalog  does  not  extend  to  periods  below  45  sec,  we  apply 
a  cos2  taper  to  the  spectra  of  observed  and  theoretical  seismograms  in  the  period  range  from 
45  to  60  sec,  in  order  to  avoid  the  ringing  associated  with  an  abrupt  truncation.  Because  of  the 
response  of  the  instruments  used,  this  leads  to  a  rather  narrow-band  signal  with  the  maximum 
energy  at  about  55  sec. 

The  signals  used  span  the  time  from  the  arrival  of  the  P-wave  until  the  arrival  of  the 
fundamental-mode  surface  waves:  Rayleigh,  for  the  vertical  and  longitudinal  components;  Love, 
for  the  transverse  component.  An  average  Earth  model  would  not  be  useful  in  the  analysis  of 
surface  waves  as,  in  the  period  range  considered  here,  their  dispersion  and  amplitudes  are 
too  much  perturbed  by  lateral  heterogeneity. 

The  segments  of  records  with  obvious  glitches  are  eliminated,  and  a  preliminary  inversion 
is  performed  to  examine  the  data  for  consistency;  this  allows  us  to  discover,  for  example,  re¬ 
versed  polarities  of  components  and  erroneous  multiplexing.  After  these  corrections  are  made, 
the  data  are  inverted  again  to  obtain  the  starting  solution:  M.°, 

At  this  stage  the  inverse  problem  is  expanded  to  allow  for  perturbations  in  the  source  radius 
r  ,  source  co-latitude  0  ,  and  longitude  X  ,  as  well  as  the  origin  time  tn.  For  the  ith  iteration, 
the  equations  of  condition  associated  with  the  k1  recording  are: 

uk(«  -  •  (t)  Mj0-11  =  6Mt  +  bk(t)  6rs 

+  ck(t)  69s  +  dR(t)  6\s  +  ek(t)  atQ  ;  (111-24) 
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with  the  implied  summation  over  i  (i  =1,2 


6),  and  where 


bk(t)  =  «/i)rlA1|{'1)(t)l  M,'-1"11 

s 


ck(t)  =  a/ao  (A1^"1)(t)]o=0  . 

s 

dk(t)  =  sinOs  a/ax  ; 

ek(t)  =  a/8t  lAk|-1)(t)l  . 

Coefficients  needed  to  evaluate  the  kernel  functions  b,  c,  d,  and  e  are  given  in  Sec.  F  above. 

Equivalent  expressions  could  also  be  derived  by  the  appropriate  differ  entiation  of  Eqs.  (2.2.30) 

1 3 

and  (2.2.31)  of  Gilbert  and  Dziewonski.  The  inverse  problem  is  usually  solved  in  the  frequency 
domain. 

Our  experience  shows  that  the  process  converges  rapidly;  almost  without  exception,  no 
substantial  decrease  in  the  rms  error  or  change  in  the  solution  is  encountered  after  third  itera¬ 
tion,  even  if  the  horizontal  shift  of  the  source  was  as  large  as  140  km,  encountered  in  me  case 
of  the  Bouvet  Island  earthquake  of  18  December  1978. 

Numerical  results  of  our  analysis  of  eight  earthquakes  ranging  in  depth  from  10  to  600  km 
and  in  seismic  moment  from  5  X  1024  to  4  x  102^  dyn-cm  are  listed  in  Table  III -11.  We  com¬ 
pare  the  locations  obtained  from  the  standard  methods.  Typically,  data  from  12  to  14  stations 
were  used  in  the  analysis. 

The  resulting  moment  tensors  are  represented  here  in  terms  of  the  moments  and  directions 
of  the  principal  axes  (tension  is  positive,  compression  is  negative).  The  solutions  have  been 
constrained  to  have  a  vanishing  trace  of  the  tensor  (Mrr  +  Mqq  +  =  0);  in  our  experiments 

we  have  found  no  evidence  that  the  isotropic  component  is  necessary  to  satisfy  the  data.  This 
does  not  mean  that  this  component  is  not  present;  our  data  are  so  clearly  dominated  by  the  shear 
energy  that,  at  least  in  this  stage  of  the  experiment,  attempts  to  recover  the  monopole  com¬ 
ponent  of  the  source  do  not  seem  warranted. 

It  was  not  required,  however,  that  the  source  should  be  of  the  double-couple  type.  Yet,  in 
six  out  of  eight  solutions  and  for  all  shallow  sources,  one  of  the  eigenvalues  is,  for  practical 
purposes,  zero.  This  means  that  the  data  demand  that  the  source  be  a  double  couple.  It  is  also 
of  interest  to  note  that  the  intermediate  eigenvalue  ("null"  axis,  for  a  double -couple  solution) 
generally  decreased  as  a  result  of  iteration  for  the  hypocentral  location.  This  may  mean  that 
bias  in  moment-tensor  solutions  could  be  introduced  by  using  the  standard  locations  or  by  lateral 
heterogeneities. 

In  the  remaining  two  cases,  for  which  the  intermediate  eigenvalue  is  nonzero,  there  is  some 
evidence  that  measurable  improvement  in  fitting  synthetics  to  the  data  can  be  accomplished  by 
introducing  nonzero  stress  along  all  principal  axes  of  the  moment  tensor.  Clearly,  many  more 
deep  and  inte  mediate  events  should  be  analyzed  in  order  to  determine  whether  such  stress  pat¬ 
terns  are  consistent  and  statistically  significant.  It  is  an  issue  of  substantial  importance  to  our 
understanding  of  the  nature  of  earthquakes  at  depths  at  which  brittle  fracture  is  not  likely. 

Figure  III-10  shows  the  distribution  of  stations  used  to  analyze  a  600-km-deep  shock  under 
the  Fiji  Islands.  As  for  most  of  the  earthquakes  in  the  circum-Pacific  belt,  the  azimuthal  dis¬ 
tribution  of  stations  is  not  well  balanced;  9  out  of  14  stations  are  in  a  single  azimuthal  quadrant. 
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Addition  of  tin-  planned  digital  WWSSN  stations  would  significantly  contribute  to  an  .mproved 
azimuthal  coverage.  The  hypocentral  parameters  of  the  event  of  2*1  April  1979  changed  very 
little.  roughly  5  km  in  location  and  ■*  i.2  sec  in  origin  time.  In  Figs.  III-ll(a)  and  (b)  we  com¬ 
part'  the  observed  seismograms  with  the  synthetics.  In  most  eases,  the  fit  is  very  good.  We 
wish  to  draw  attention  to  the  traces  of  the  vertical  component  of  station  KONO;  the  analyzed 
part  of  the  recording,  delineated  by  broken  vertical  lines,  is  45  min.  long,  and  the  synthetics 
reproduce  the  observed  record  in  all  important  details.  The  same  is  true  in  respect  to  nearly 
all  vertical  and  transverse  components,  agreement  for  SNZO,  NWAO,  and  ANMO  (transverse 
component)  is  spectacular.  However,  the  synthetics  for  longitudinal  components  predict  wave- 
trains  that,  at  late  times,  have  systematically  lower  amplitudes.  This  cannot  be  explained  by 
an  error  in  source  mechanism  or  the  Q-strueture  of  the  model,  since  we  would  not  then  expect 
a  satisfactory  agreement  for  vertical  and  transverse  components.  We  think  that  the  explanation 
must  be  related  to  differences  between  the  elastic  parameters  of  the  upper  mantle  in  the  real 
Far th  and  in  the  model  1066B.  The  effect  is  clearly  a  noticeable  one;  it  is  possible  to  think  that 
studies  such  as  this  one  will  lead  to  the  refinement  of  the  Earth  structure  in  addition  to  con¬ 
tributing  to  our  knowledge  of  global  seismicity. 

Figures  III-12(a)  and  (b)  shows  selected  seismograms  and  synthetics  for  an  intermediate- 
depth  shock  of  23  September  1978.  The  epicenter  is  not  too  distant  from  that  discussed  before 
and,  consequently,  station  coverage  is  very  similar  to  that  shown  in  Fig.  III-10.  The  change 
in  epicentral  coordinates  is  much  greater  in  this  case;  depth  decreased  by  8  km  and,  although 
the  latitude  remained  the  same,  longitude  decreased  by  0.46’;  also,  the  origin  time  increased 
by  9.3  sec.  If  this  difference  is  not  an  artifact  of  lateral  heterogeneities  (the  fact  that  the  loca¬ 
tion  of  the  24  April  1979  event  remained  unchanged  would  tend  to  indicate  that  this  is  not  the 
case)  or  poor  standard  location,  then  one  would  have  to  infer  that  this  earthquake  of  seismic 
moment  1.6  X  10 ^  dyn-cm  had  substantial  source  dimensions  (~100  km)  and  a  duration  of 
roughly  20  sec. 

Comparison  of  observed  seismograms  with  the  synthetics  gives  an  opportunity  to  appreciate 
path  dependence  of  the  quality  of  fit.  Most  of  the  stations  fit  well;  an  important  exception  is 
represented  by  stations  MAIO  and  KAAO;  the  paths  are  very  similar  for  these  two  stations 
and  they  cross  the  most  tectonically  disturbed  region  of  Asia.  Although  the  observed  and  syn¬ 
thetic  seismograms  show  good  agreement  for  early  phases  (traveling  through  the  lower  mantle), 
the  multiply  reflected  Sy  waves  show  substantial  differences;  about  35  min.  after  the  origin  time 
the  synthetics  are  180°  out  of  phase  with  respect  to  observations  for  the  vertical  and  longitudinal 
components.  It  is  interesting  to  note  that  the  effect  on  waves  is  much  less;  here  the  agree¬ 
ment  is  rather  good  throughout.  Observations  at  CHTO,  closer  to  the  epicenter  but  only  10° 
different  in  azimuth  at  the  source,  are  in  excellent  agreement  with  the  synthetics,  confirming 
the  evidence  for  anomalous  properties  of  the  remainder  of  the  path. 

Each  of  the  eight  analyses  brought  to  light  some  new  and  unexpected,  at  least  by  us,  infor¬ 
mation.  Because  of  the  format  of  this  report,  we  cannot  provide  a  detailed  description  of  the 
results.  We  feel  that  the  approach  outlined  here  can  be  efficiently  used  to  learn  more  about 
the  earthquakes  and  the  Earth.  Teh- An  Chou* 

A.  M.  Dziewonskj 


*  Department  of  Geological  Sciences,  Harvard  University. 
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Fig.  III-l.  Short-period  waveforms  observed  at  four  SRO/ASRO  stations 
from  events  7,  8,  and  9  (see  Table  III- 1 ). 
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Fig.IIl-2.  Rayleigh-wave  amplitudes  as  a  function  of  azimuth: 
(a)  peak  amplitude  on  seismogram;  (b)  30 -sec  period  spectral 
amplitudes,  corrceted  for  geometrical  spreading  and  attenuation 
of  0.75  x  lOy  km;  and  (c)  same  as  (b)  but  20-sec  period. 
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Kig.  Ill- 3 .  Rayleigh-wave  seismograms 
for  events  7,  8,  and  9  at  KAAO  (Kabul) 
and  MAJO  (Matsushiro)  showing  phase 
reversal  of  seismograms  for  event  8. 
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Fig.  III-4,  Peak  Love-wave  amplitudes  as  a  function  of  azimuth, 
measured  from  seismograms. 
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Fig.  III-7.  Values  of  Al,  and  Bl  at  a  variety  of  source  depths.  Al  indicated 
by  solid  line,  Bj^  by  dashed  line. 
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Fig.  Ill -8.  Range  of  eigenvalues  obtained  as  a  function  of  source 
depth  at  which  synthetic  data  are  generated. 
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Fig.  III-9.  Variation  of  (a)  Rayleigh-wave  and  (b)  Love-wave 
spectral  amplitudes  with  azimuth,  as  predicted  from  moment- 
tensor  solution  obtained  for  events  7,  8,  and  9. 
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Fie.  III-l  1  (a-b) .  Deep  (600-km)  earthquake  of  M  April  *9 79:  comparison 
of  observed  and  y^^^ommonto^oth  observed  and  synthetic 

^trTerca'dastd"  Ss  de-i^ta^^w^ 

lyzed.  Theoretical  seismograms  were  not  computed  o  m  (1)  2-1  * 

val25  aPrinCiPaZ3V/ai6U8S  mo  4  xTo^d^cm,  5/76;  (!)  -2.5  X  1025  dyn-cm 
66/334!"  tension^ is  ^^sitive'.^  Notice  gTitch  associated  with  P-wave  arnval 
on  record  for  station  ANMO~Z. 
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Fig.  Ill —1 Z (a -b) .  Same  as  Fig.  Ill -1 1  but  for  event  of  23  September  1978.  Principal 
values  and  directions  of  axes  (plunge/Az)  are:  (1)  1.4  X  102°  dyn-cm,  79/154; 
(2)  0.0  dyn-cm,  11/338;  (3)  -1.4  X  1026  dyn-cm,  0/248;  tension  is  positive. 
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IV.  COMPUTER  SYSTEMS 


A.  UNIX  SYSTEM  ENHANCEMENTS 

Several  efforts  are  being  conducted  to  improve  the  computational  support  available  to  UNIX 
users  at  Lincoln  Laboratory.  The  major  projects  are  the  installation  of  a  new  release  of  Bell 
Laboratories'  UNIX  operating  system;  coordination  of  software  maintenance  functions  at  existing 
UNIX  groups  within  the  Laboratory;  and  improvement  of  ARPANET  software. 

1.  Installation  of  Seventh- Edition  UNIX  (V7) 

We  are  presently  completing  installation  of  the  seventh-edition  release  of  the  UNIX  operat¬ 
ing  system  (V7)  on  four  computer  systems:  a  PDP-ll/45  and  a  PDP-ll/50  used  by  Group  22 
largely  for  computational  support  seismology  research,  a  PDP-ll/45  used  by  Group  24  for  speech 
research,  and  a  PDP-ll/70  operated  jointly  by  the  two  Groups.  The  V7  provides  a  number  of 
features  which  were  unavailable  in  previous  releases.  It  is  available  for  both  the  PDP-11  and  the 
VAX  CPU  families,  and  permits  substantially  larger  data  sets  than  the  previous  versions.  In 
addition,  a  number  of  bugs  have  been  fixed  and  overall  performance  has  been  improved. 

The  user  environment  of  V7  is  also  significantly  improved.  Bell  provides  a  dialect  of  For¬ 
tran  77  which  allows  convenient  mixing  and  interfacing  of  Fortran  and  C  routines.  A  number  of 
programs  for  performing  routine  data-processing  tasks  have  been  added  and  the  text-processing 
package  has  been  enhanced. 

To  speed  the  introduction  of  V7  programs  and  to  smooth  the  transition  to  the  new  system, 
we  made  a  few  modifications  to  our  existing  sixth-edition  (V6)  kernel.  The  modified  kernel  pro¬ 
vides  the  additional  services  required  to  execute  user  programs  written  for  either  V6  or  V7. 

This  has  allowed  us  to  run  all  of  the  V7  programs  and  still  maintain  access  to  the  ARPANET, 
for  which  no  V7  compatible  software  was  available.  Also,  it  has  provided  a  period  during  which 
old  programs  continue  to  function  and  may  be  modified  to  run  on  a  standard  V7  system.  As  soon 
as  final  debugging  of  the  V7  ARPANET  software  is  finished,  we  will  install  the  new  kernel,  which 
provides  large  file  systems  and  other  valuable  although  noncritical  features. 

2.  Site  Compatibility 

The  modified  V6  kernel  has  been  installed  on  the  four  existing  UNIX  systems  at  Lincoln. 

This  has  eliminated  the  growing  divergence  between  the  software  at  each  site  and  has  reduced 
the  amount  of  duplicated  effort  expended  in  maintaining  the  systems.  Since  all  software  will  run 
on  any  site  without  modification,  corrections  and  enhancements  done  at  any  site  may  be  quickly 
distributed  to  the  other  sites.  This  compatibility  has  been  enhanced  by  the  use  of  tables  which 
define  site-specific  parameters  (e.  g.,  ARPANET  host  names).  These  parameters  are  read  when 
programs  are  executed  rather  than  being  fixed  at  compile  time.  This  substantially  reduces  the 
effort  involved  in  distributing  a  program,  while  introducing  no  noticeable  change  in  the  program's 
performance. 

3.  New  ARPANET  Software 

As  the  final  phase  of  the  conversion,  the  BBN  ARPANET  software  is  being  modified  to  inter¬ 
face  with  the  V7  kernel.  This  will  correct  many  bugs  in  the  earlier  code  obtained  from  the  Uni¬ 
versity  of  Illinois,  as  well  as  introducing  the  new  9 6 -bit  leader  formats.  This  effort  is  in  the 
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final  debugging  stages.  When  the  installation  is  completed,  the  V7  kernel  with  network  interface 
will  replace  the  present  modified  \<>  kernel  at  all  sites. 

K.  J.  Schroder 

R.  Graphpac 

Graphpac  is  a  package  of  over  40  subroutines  and  functions  which  are  used  to  produce  a 
graphics  display  on  the  Tektronix  4014  terminal.  The  routines  can  be  called  from  both  C  and 
Fortran  programs.  Lower-level  Graphpac  routines  include  drawing  a  line  from  one  point  on  the 
screen  to  another,  moving  the  cursor  from  one  point  to  another,  outputting  text,  and  erasing  the 
screen.  Higher-level  routines  include  drawing  linear-  and  log-scaled  axes,  plotting  arrays  of 
data,  drawing  polygons,  and  determining  the  position  of  the  cross  hairs.  Graphpac  was  written 
to  improve  the  speed  of  generating  graphics,  and  to  provide  higher-level  graphics  functions. 

By  not  using  a  daemon,  and  sending  output  immediately  to  the  user's  terminal,  speed  w*as  gained 
over  existing  graphics  routines.  The  higher-level  Graphpac  functions  have  proved  to  be  useful 
in  several  applications  programs  requiring  axes,  plotting  of  arrays,  and  other  higher-level 
functions. 

An  important  part  of  Graphpac  is  the  concept  of  windows.  The  "uv  window,"  set  by  the  user, 
is  a  rectangular  subset  of  the  Tektronix  screen.  Once  the  user  sets  the  uv  window,  he  can  draw 
and  move  the  cursor  in  screen  coordinates.  A  full-screen  uv  window  has  abscissas  from  0  to 
4095,  and  ordinates  from  0  to  3119.  However,  the  user's  data  values  usually  do  not  correspond 
with  the  values  of  the  screen  coordinates.  The  "xy  window,"  set  by  the  user,  is  the  rectangular 
subset  of  the  Cartesian  coordinate  plane  where  his  data  resides.  Once  the  xy  window  and  the  uv 
window  have  been  set,  Graphpac  will  automatically  map  coordinates  in  the  xy  window  to  their 
corresponding  values  in  the  uv  window.  The  user  can  specify  whether  the  abscissas  and  ordi¬ 
nates  in  the  xy  window  are  linear  or  log  scaled. 

Graphpac  also  permits  clipping  to  be  performed.  When  clipping  i.’  activated,  subsequent 
commands  to  draw  lines  will  display  only  the  portions  of  the  lines  wh'.jh  reside  within  the  uv 
window.  Portions  of  lines  which  are  outside  the  uv  window  are  "clipped"  and  not  displayed. 
Coordinates  in  the  xy  window  are  first  mapped  to  the  uv  window  before  clipping  is  done. 

P.  T.  Crames 


C.  CONVERSION  TO  FORTRAN  77 

Our  Group  has  been  bilingual  with  respect  to  programming  languages  since  we  began  using 
UNIX.  In  general,  Fortran  is  preferred  by  the  seismologists  because  they  are  already  familiar 
with  it,  have  existing  Fortran  programs  written  for  use  on  other  systems,  and  would  like  to 
write  in  a  language  that  is  portable  to  other  machines  with  a  minimum  number  of  changes.  Most 
of  the  Group's  programmers,  on  the  other  hand,  prefer  the  standard  UNIX  language  "C"  because 
of  its  structured  philosophy,  ease  of  interactive  use,  flexibility,  efficiency,  and  code  compact¬ 
ness.  The  use  of  two  languages  for  applications  programs  has  resulted  in  frequent  unnecessary 
duplication  of  existing  software  and  necessitated  the  development  of  two  versions  of  several  sub¬ 
routine  packages,  with  the  attendant  difficulties  of  parallel  maintenance. 

Se  venth-edition  UNIX  is  supplied  with  a  Fortran  compiler  that  closely  matches  the  new  For¬ 
tran  standard,  commonly  known  as  Fortran  77.  This  new  compiler  has  the  attractive  feature  of 
generating  code  that  is  compatible  with  the  code  generated  from  C.  Interface  subroutines  were 
written  in  C  for  several  libraries  of  subprograms  having  seismic  applications,  and  a  conversion 
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document  was  distributed  which  suggested  equivalent  alternatives  to  those  features  of  Princeton 
UNIX  Fortran  which  wc  had  been  using  that  would  not  work  with  the  new  compiler. 

‘  After  considerable  use.  we  have  found  that  the  new  compiler  is  slower  and  produces  muc 
longer  executable  code  because  of  a  rather  cumbersome  implementation  of  formatted  I/O;  how¬ 
ever  the  advantages  of  C  language  compatibility,  greater  availa.le  program  data  space, 
diagnostics,  and  a  wider  range  of  intrinsic  functions  make  the  use  of  the  Fortran  77  compiler 


preferable  for  many  applications. 

K  D.  A.  Bach 
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GLOSSARY 


AA 

ARPANET 

ASRO 

Automatic  Association 

DARPA  Computer  Network 

Upgraded  HGLP  Station 

BBN 

Bolt,  Beranek  and  Newman,  Inc. 

CLVD 

CPU 

Compensated  Linear  Vector  Dipole 

Control  and  Processing  Unit 

DADS 

DARPA 

DEC 

Data  Analysis  and  Display  System 

Defense  Advanced  Research  Projects  Agency 
Digital  Equipment  Corporation 

HGLP 

High -Gain,  Long-Period  (Station) 

ISC 

ISDE 

ISM 

International  Seismological  Ccf'* 

International  Seismic  Data  Exchange 

International  Seismic  Month 

LRU 

Least  Recently  Used 

NMRO 

Nuclear  Monitoring  Research  Office 

PDE 

Preliminary  Determination  of  Epicenter 

SAS 

SATS 

SDC 

SRO 

STAP 

Seismic  Analysis  Station 

Semiannual  Technical  Summary 

Seismic  Data  Center 

Seismic  Research  Observatory 

Short-Term  Average  Power 

UNIX 

USGS 

(Trade-Mark)  Bell  Laboratories  Operating  System 
U.S.  Geological  Survey 

VSC 

Vela  Seismological  Center 

WMO 

WWSSN 

World  Meteorological  Organization 

World-Wide  Standard  Seismograph  Network 
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